Calculation of repertoire metrics
VDJTOOLS的CalcDiversityStats函数提供了直接计算功能，详情参考免疫组库分析实战/mixcr+vdjtools+R实现，使用vdjtools进行免疫组库分析.
本文以TRUST4的结果为例（基本都大同小异），展示如何计算免疫组库的Richness/Diversity/Clonality。
trust4.report
-
Diversity：香农-威纳指数（Shannon-WienerIndex ，H） 可以回答样本各物种均一性如何，一般选择以e为底的对数计算。
在香农-威纳指数中，包含着两个成分：①种数（Richness）；②各种间个体分配的均匀性（Evenness）
Shannon-WienerIndex Richness : 每个样本中的Number of Unique clone。
Evenness：各种之间，个体分配越均匀，H值就越大。如果每一个体都属于不同的种，多样性指数就最大；如果每一个体都属于同一种，则其多样性指数就最小。可以通过估计群落的理论上的最大多样性指数（Hmax），然后以实际的多样性指数对Hmax的比率，从而获得均匀性指数，即为标准化的香农指数。
Evenness=H/Hmax，当每一个Clone出现的概率相同时即可得到Hmax，此时，P=1/N，Hmax=-1/N*log(1/N)*N=-log(1/N)=log(N), Evenness=Diversity/log(Richness) 。Clonality: Clonality = 1 – Evenness = 1 - Diversity/log(Richness)。
参考：Shannon–Wiener Index
代码实现（R）
Shannon_index <- function(list_p){
sum=0
for (p in list_p){
sum <- p*log(p) + sum
}
H <- -sum
return(H)
}
clonality <- function(diversity,richness){
E<- diversity/log(richness)
C <- 1-E
return(C)
}
sample_id <- c(paste0("sample",1:10))
results <- data.frame(sample_id)
for (i in 1:length(sample_id)){
df <- read.csv(file = paste0(sample_id[i],".txt"),
stringsAsFactors = F,sep = "\t")
list_p <- df[,"frequency"]
results[i,"richness"] <- nrow(df)
results[i,"diversity"] <- Shannon_index(list_p)
results[i,"clonality"] <- clonality(Shannon_index(list_p),nrow(df))
}
write.table(results,filename,sep = "\t",row.names = F,quote = F)
