---
title: "Running LACE"
author: 
  - Daniele Ramazzotti
  - Fabrizio Angaroni
  - Davide Maspero
  - Alex Graudenxi
  - Luca De Sano
  - Gianluca Ascolani
date: "`r format(Sys.time(), '%B %d, %Y')`"
graphics: yes
package: LACE
output:
  BiocStyle::html_document:
    toc_float: true
vignette: >
  %\VignetteIndexEntry{Running LACE}
  %\VignetteEngine{knitr::rmarkdown}
  %\VignetteEncoding{UTF-8}
  %\VignetteDepends{LACE,BiocStyle}
---


## Setup the environement

We now present an example of longitudinal analysis of cancer evolution with *LACE* using single-cell data obtained from Rambow, Florian, et al. "Toward minimal residual disease-directed therapy in melanoma." Cell 174.4 (2018): 843-855. 

As a first step, we load single-cell data for a skin melanoma. The data comprises point mutations for four 
time points: (1) before treatment, (2) 4 days treatment, (3) 28 days treatment and finally (4) 57 days treatment. 

```{r}
library("LACE")
data(longitudinal_sc_variants)
names(longitudinal_sc_variants)
```

We setup the main parameter in oder to perform the inference. First of all, as the four data points may potentially provide sequencing for an unbalanced 
number of cells, we weight each time point proportionally to the sample sizes as follow. We refer to the paper for details. 

```{r}
lik_weights = c(0.2308772,0.2554386,0.2701754,0.2435088)
```

The second main parameter to be defined as input is represented by the false positive and false negative error rates, i.e., alpha and beta. We can specify a 
different rate per time point as a list of rates. When multiple set of rates are provided, *LACE* performs a grid search in order to estimate the best set of error rates. 

```{r}
alpha = list()
alpha[[1]] = c(0.02,0.01,0.01,0.01)
alpha[[2]] = c(0.10,0.05,0.05,0.05)
beta = list()
beta[[1]] = c(0.10,0.05,0.05,0.05)
beta[[2]] = c(0.10,0.05,0.05,0.05)
head(alpha)
head(beta)
```

## Inference

We can now perform the inference as follow. 

```{r}
inference = LACE(D = longitudinal_sc_variants, 
    lik_w = lik_weights, 
    alpha = alpha, 
    beta = beta, 
    keep_equivalent = TRUE, 
    num_rs = 5, 
    num_iter = 10, 
    n_try_bs = 5, 
    num_processes = NA, 
    seed = 12345, 
    verbose = FALSE)
```

We notice that the inference resulting on the command above should be considered only as an example; the parameters num_rs, num_iter and n_try_bs representing the number of 
steps perfomed during the inference are downscaled to reduce execution time. We refer to the Manual for discussion on default values. We provide within the package results 
of inferences performed with correct parameters as RData. 

```{r}
data(inference)
print(names(inference))
```

*LACE* returns a list of nine elements as results. Namely, B and C provide respectively the maximum likelihood longitudinal tree and cells attachments; corrected_genotypes the corrected genotypes, clones_prevalence, the estimated prevalence of any observed clone; relative_likelihoods and joint_likelihood the estimated likelihoods for each time point and the weighted likelihood; clones_summary provide a summary of association of mutations to clones. In equivalent_solutions, solutions (B and C) with likelihood equivalent to the best solution are returned; notice that in the example we disabled this feature by 
setting equivalent_solutions parameter to FALSE. Finally, error rates provide the best error rates (alpha and beta) as estimated by the grid search. 


## Plot

We can plot the inferred model using the function longitudinal.tree.plot. 

```{r fig.width=12, fig.height=8, warning=FALSE, fig.cap="Inferred model"}
clone_labels = c("ARPC2","PRAME","HNRNPC","COL1A2","RPL5","CCT8")
longitudinal.tree = longitudinal.tree.plot(inference = inference, 
                                           labels_show = "clones", 
                                           clone_labels = clone_labels, 
                                           legend_position = "topright")
```
