TCR_analysis.Rmd

---
title: "TCR clonotype analysis"
author: "Benoit P Nicolet"
date: "16-3-21"
output: html_document
---

```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = TRUE)

#BiocManager::install("dtplyr")
#install.packages("dtplyr") # skip this if you already installed devtools
#devtools::install_github("immunomind/immunarch")

#BiocManager::install("immunarch")

library(plyr)
library(dplyr)
library(dtplyr)
library(stringr)
library(ggplot2)
library(immunarch)

knitr::opts_knit$set("/path-to-analysis-folder/")

```


```{r importing data}

#data(immdata) # This will use the data of the package. 

immdata <- repLoad("/path-to-mixcr-output/")  # Replace it with the path to your data. Immunarch automatically detects the file format. 
# <3 lovely function.... 

head(immdata) # list of dataframe

```


```{r visualizing stuff}

## Clonotype versus CDR3 region length ##
repExplore(immdata$data, "lens") %>% vis() #+ theme(legend.text = element_blank())  # Visualise the length distribution of CDR3

## Colonotype abundance group per sample ##
repClonality(immdata$data, "homeo") %>% vis()  # Visualise the relative abundance of clonotypes
# Keep in mind here that the sample sequencing depth will greatly influence the clonotype abundance...!


## Shared clonotypes between samples ##
repOverlap(immdata$data) %>% vis()  # Build the heatmap of public clonotypes shared between repertoires


## TCR V genes used in sample 1 ##
geneUsage(immdata$data[[1]]) %>% vis()  # Visualise the V-gene distribution for the first repertoire  (change number for the sample number)

## diversity analysis ##
repDiversity(immdata$data) %>% vis(.meta = immdata$meta)  # Visualise the Chao1 diversity of repertoires, grouped by the patient status
# here it uses Chao1, but you can change the method using .method= ......

## Waterfall plot of colotype and abundance ## 
repExplore(immdata$data, .method = "count") %>% vis() ## pretty pretty : )

## Rarefaction analysis ## 
repDiversity(immdata$data, "raref", .verbose = F) %>% vis() # computationally intense for sample size >6 ...

## tracking shared clonotype throughout samples ##
trackClonotypes(immdata$data,list(2, 10), .col = "nt") %>% vis() + theme(legend.text = element_blank())
# here in list(x, y), x is the sample number and y is the number of top TCR sequences to consider. 

```