---
title: "Famdenovo.brca Example"
output:
  html_document:
    df_print: kable
---

```{r setup, include=FALSE}
knitr::opts_chunk$set(fig.width=12, fig.height=8) 
#note, Famdenovo.BRCA requires BayesMendel, which can be obtained here: https://projects.iq.harvard.edu/bayesmendel/bayesmendel-r-package
require(BayesMendel)
source("modified.BRCAPRO.combined.R")
source("runPeeling.combined.R")
source("Famdenovo.BRCA.combined.R")
library(ggplot2)
library(dplyr)
library(kinship2)

```

## Pedigree Format


+-----------------------+-----------------------------------+
| Column Name           |Content|
+=======================+===================================+
| FamilyID              |Family identifier|
+-----------------------+-----------------------------------+
| BRCA1carrier          |BRCA1 carrier status (NA=unknown, 0=confirmed normal, 1=confirmed carrier)|
+-----------------------+-----------------------------------+
| BRCA2carrier          |BRCA2 carrier status ( NA=unknown,0=confirmed normal, 1=confirmed carrier)|
+-----------------------+-----------------------------------+
| ID                    |Member identifier|
+-----------------------+-----------------------------------+
| Gender                |Gender (0=female, 1=male)|
+-----------------------+-----------------------------------+
| FatherID              |Father’s identifier number|
+-----------------------+-----------------------------------+
|   MotherID            |Mother’s identifier number|
+-----------------------+-----------------------------------+
|   AffectedBreast      |Breast cancer status (0=no cancer, 1=breast cancer, one breast involved; 2=bilateral breast cancer, NA=unknown status)|
+-----------------------+-----------------------------------+
|   AffectedOvary       |Ovarian cancer status (0=no cancer, 1=ovarian cancer, NA=unknown status)|
+-----------------------+-----------------------------------+
|   AgeBreast           |Age of onset of breast cancer if a breast cancer case.Current age or age of death if not a breast cancer case.NA if there is no age information|
+-----------------------+-----------------------------------+
|  AgeOvary             |Age of onset of ovarian cancer if an ovarian cancer case. Current age or age of death if not an ovarian cancer case.NA if there is no age information|
+-----------------------+-----------------------------------+
|AgeBreastContralateral |Age at onset of breast cancer, second breast.Only for members with breast cancer status=2. For the rest enter a 0.|
+-----------------------+-----------------------------------+
| Twins                  Identifies siblings who are identical twins. Enter 0 for everyone. Functionality to handle twins will be developed in future.|
+-----------------------+-----------------------------------+
| ethnic                |Identifies the ethnicity of each family member.Enter “nonAJ”, “AJ”, “Italian”, “Other” or NA (as recognized by is.na() function).|
+-----------------------+-----------------------------------+

## Example Pedigree Data

```{r,echo=F}
sample_fam<-data.frame(FamilyID=c(rep(1,21)),BRCA1carrier=c(rep(0,21)), BRCA2carrier=c(1,rep(NA,20)), ID=c(1:21), Gender=c(0,0,1,0,0,1,1,0,0,0,1,0,1,0,0,0,0,0,1,1,1), FatherID=c(3,11,13,21,21,21,3,3,3,0,0,0,0,11,11,11,13,13,13,13,0), MotherID=c(2,10,12,1,1,1,2,2,2,0,0,0,0,10,10,10,12,12,12,12,0), AffectedBreast=c(1,1,0,0,0,0,0,0,0,0,0,0,0,0,2,0,0,0,0,0,0), AffectedOvary=c(0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0), AgeBreast=c(43,28,67,15,18,13,40,41,37,77,55,93,82,9,44,53,70,61,73,55,47),AgeOvary=c(44,53,67,15,18,13,40,41,37,77,55,93,82,9,44,53,70,61,73,55,47), AgeBreastContralateral=c(0,0,0,0,0,0,0,0,0,0,0,0,0,0,44,0,0,0,0,0,0), Twins=c(rep(0,21)), ethnic=rep("nonAJ",21))


sample_fam



```


```{r,echo=F}

  sample_fam<-sample_fam%>%rowwise()%>%mutate(leastage=min(AgeBreast,AgeOvary))%>%ungroup()
  IDdf<- data.frame("ID"=sample_fam$ID,"AgeID"= paste(sample_fam$ID, sample_fam$leastage, sep = "\n Age:"),stringsAsFactors = FALSE)
  for (i in 1:length(IDdf$ID)){
    sample_fam$ID[i]<-IDdf$AgeID[i]
    if(sample_fam$FatherID[i]!=0){
    fid<-IDdf%>%filter(ID==sample_fam$FatherID[i])
    sample_fam$FatherID[i]<-as.character(fid$AgeID)}
    if(sample_fam$MotherID[i]!=0){
    mid<-IDdf%>%filter(ID==sample_fam$MotherID[i])
    sample_fam$MotherID[i]<-as.character(mid$AgeID)}
  }



sample_fam$Gender[which(sample_fam$Gender==0)]<-2
sample_fam$AffectedBreast[which(sample_fam$AffectedBreast==2)]<-1
sample_fam$BRCA2carrier[is.na(sample_fam$BRCA2carrier)]=0  
aff2 <- data.frame(BreastCancer=sample_fam$AffectedBreast, BRCA =sample_fam$BRCA2carrier)
 
ped<-pedigree(id=sample_fam$ID,dadid=sample_fam$FatherID,momid=sample_fam$MotherID,sex = sample_fam$Gender,affected = as.matrix(aff2),missid = 0)
plot.pedigree(ped)
pedigree.legend(ped, location = "bottomright", radius = 0.1)


```

## Calling Famdenovo.BRCA


Famdenovo.BRCA() takes three inputs:

* Pedigree Data, as described above
 * _De_ _novo_ ratio for _BRCA1_
 * _De_ _novo_ ratio for _BRCA2_

```{r, include=FALSE}
sample_fam<-data.frame(FamilyID=c(rep(1,21)),BRCA1carrier=c(rep(0,21)), BRCA2carrier=c(1,rep(NA,20)), ID=c(1:21), Gender=c(0,0,1,0,0,1,1,0,0,0,1,0,1,0,0,0,0,0,1,1,1), FatherID=c(3,11,13,21,21,21,3,3,3,0,0,0,0,11,11,11,13,13,13,13,0), MotherID=c(2,10,12,1,1,1,2,2,2,0,0,0,0,10,10,10,12,12,12,12,0), AffectedBreast=c(1,1,0,0,0,0,0,0,0,0,0,0,0,0,2,0,0,0,0,0,0), AffectedOvary=c(0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0), AgeBreast=c(43,28,67,15,18,13,40,41,37,77,55,93,82,9,44,53,70,61,73,55,47),AgeOvary=c(44,53,67,15,18,13,40,41,37,77,55,93,82,9,44,53,70,61,73,55,47), AgeBreastContralateral=c(0,0,0,0,0,0,0,0,0,0,0,0,0,0,44,0,0,0,0,0,0), Twins=c(rep(0,21)), ethnic=rep("nonAJ",21))

```
```{r}
Famdenovo.BRCA(sample_fam,denovo_ratio1 = .15, denovo_ratio2 = .15)



```

For this sample family, the one carrier, ID = 1, has a _BRCA2_ mutation. This mutation is not predicted to be _de_ _novo_, with a probability of 0.0029








## _De novo_ Example
Here we make a few changes to the pedigree, to make it a _de_ _novo_ example. We do this by making the earlier generations less likely to carry mutation, and later generations more likely. We also change ID 4 to a _BRCA2_ carrier as an additional example of a familial mutation.  

```{r,echo=F}
sample_fam<-data.frame(FamilyID=c(rep(1,21)),BRCA1carrier=c(rep(0,21)), BRCA2carrier=c(1,0,0,1,rep(NA,17)), ID=c(1:21), Gender=c(0,0,1,0,0,1,1,0,0,0,1,0,1,0,0,0,0,0,1,1,1), FatherID=c(3,11,13,21,21,21,3,3,3,0,0,0,0,11,11,11,13,13,13,13,0), MotherID=c(2,10,12,1,1,1,2,2,2,0,0,0,0,10,10,10,12,12,12,12,0), AffectedBreast=c(1,0,0,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0), AffectedOvary=c(0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0), AgeBreast=c(43,80,70,25,35,13,40,41,37,77,75,93,82,9,54,53,70,61,73,55,47),AgeOvary=c(44,80,70,30,46,13,40,41,37,77,75,93,82,9,65,53,70,61,73,55,47), AgeBreastContralateral=c(0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0), Twins=c(rep(0,21)), ethnic=rep("nonAJ",21))




  sample_fam<-sample_fam%>%rowwise()%>%mutate(leastage=min(AgeBreast,AgeOvary))%>%ungroup()
  IDdf<- data.frame("ID"=sample_fam$ID,"AgeID"= paste(sample_fam$ID, sample_fam$leastage, sep = "\n Age:"),stringsAsFactors = FALSE)
  for (i in 1:length(IDdf$ID)){
    sample_fam$ID[i]<-IDdf$AgeID[i]
    if(sample_fam$FatherID[i]!=0){
    fid<-IDdf%>%filter(ID==sample_fam$FatherID[i])
    sample_fam$FatherID[i]<-as.character(fid$AgeID)}
    if(sample_fam$MotherID[i]!=0){
    mid<-IDdf%>%filter(ID==sample_fam$MotherID[i])
    sample_fam$MotherID[i]<-as.character(mid$AgeID)}
  }



sample_fam$Gender[which(sample_fam$Gender==0)]<-2
sample_fam$AffectedBreast[which(sample_fam$AffectedBreast==2)]<-1
sample_fam$BRCA2carrier[is.na(sample_fam$BRCA2carrier)]=0  
aff2 <- data.frame(BreastCancer=sample_fam$AffectedBreast, BRCA =sample_fam$BRCA2carrier)
 
ped<-pedigree(id=sample_fam$ID,dadid=sample_fam$FatherID,momid=sample_fam$MotherID,sex = sample_fam$Gender,affected = as.matrix(aff2),missid = 0)
plot.pedigree(ped)
pedigree.legend(ped, location = "bottomright", radius = 0.1)



sample_fam<-data.frame(FamilyID=c(rep(1,21)),BRCA1carrier=c(rep(0,21)), BRCA2carrier=c(1,0,0,1,rep(NA,17)), ID=c(1:21), Gender=c(0,0,1,0,0,1,1,0,0,0,1,0,1,0,0,0,0,0,1,1,1), FatherID=c(3,11,13,21,21,21,3,3,3,0,0,0,0,11,11,11,13,13,13,13,0), MotherID=c(2,10,12,1,1,1,2,2,2,0,0,0,0,10,10,10,12,12,12,12,0), AffectedBreast=c(1,0,0,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0), AffectedOvary=c(0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0), AgeBreast=c(43,80,70,25,35,13,40,41,37,77,75,93,82,9,54,53,70,61,73,55,47),AgeOvary=c(44,80,70,30,46,13,40,41,37,77,75,93,82,9,65,53,70,61,73,55,47), AgeBreastContralateral=c(0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0), Twins=c(rep(0,21)), ethnic=rep("nonAJ",21))


Famdenovo.BRCA(sample_fam,denovo_ratio1 = .15, denovo_ratio2 = .15)
```



Now, both ID = 1 and 4 have _BRCA2_ mutations. For ID = 4, as we would expect, the mutation is not predicted to be _de_ _novo_ with a probability of 0.007. ID = 1 is predicted to be _de_ _novo_, with a probability of 0.510.