今天看啥  ›  专栏  ›  生信店小二

TCR-seq数据分析(二)

生信店小二  · 简书  ·  · 2021-02-19 16:25

前言

上一次《 TCR-seq数据分析一 》分享了如何初步处理TCR-seq的数据,这一次接着来说如何进一步处理数据。上一章介绍了免疫受体的V(D)J结构,当发生免疫反应时,VDJ区域的基因重排会形成很多不同的免疫受体,即免疫组库。TCR-seq数据分析的主要目的就是统计各区域基因的出现频率,即geneUsage。这个时候就轮到今天的主角上场了——immunarch是一个R包,可以用来对很多软件的TCR-seq数据如mixcr、10X等做后续的数据分析。该R包含有丰富的处理函数以及多样性的数据展示类型,用起来也是相当方便,那话不多说来看一下如何使用该R包。

代码

install.packages("immunarch")
library(immunarch)
data(immdata)
gene_stats()
##     alias                 species ighd ighj ighv igij igkj igkv iglj iglv traj
## 1      bt               BosTaurus   21    4   25    0    1    6    5   26   46
## 2      cd      CamelusDromedarius    0    0    0    0    0    0    0    0    0
## 3     clf    CanisLupusFamiliaris    0    0    0    0    0    0    0    0    0
## 4      dr              DanioRerio    7    7    0    3    0    0    0    0    0
## 5      hs             HomoSapiens   30   13  248    0    5   64    7   69   57
## 6  macmul           MacacaMulatta   24    7   19    0    4   83    5    0    0
## 7     mmc    MusMusculusCastaneus    0    0    0    0    0    4    0    0    0
## 8     mmd   MusMusculusDomesticus    0    0    0    0    0    2    0    0    0
## 9  musmus             MusMusculus   32    8  225    0    8  109    3    5   42
## 10     oa OrnithorhynchusAnatinus    3   10    0    0    0    0    0    0    0
## 11     oc    OryctolagusCuniculus   10   11   39    0    8   26    2   20    0
## 12     om      OncorhynchusMykiss    9    7    6    0    0    0    0    0    0
## 13     rn        RattusNorvegicus   30    4  113    0    6  132    2    8    0
## 14   smth   MusMusculusMolossinus    0    0    0    0    0    1    0    0    0
## 15   smth     MusMusculusMusculus    0    0    0    0    0    1    0    0    0
## 16   smth              MusSpretus    0    0    0    0    0    2    0    2    0
## 17     ss               SusScrofa    5    5   15    0    8   19    4   14    0
##    trav trbd trbj trbv trdd trdj trdv trgj trgv
## 1     0    0    0    0    5    3    0    6   15
## 2     0    0    0    0    0    0    7    2    2
## 3     0    2    8   19    0    0    0    7    8
## 4     0    0    0    0    0    0    0    0    0
## 5    60    3   14   64    3    4    6    4   10
## 6     0    2   15   58    0    0    0    0    0
## 7     0    0    0    0    0    0    0    0    0
## 8     0    0    0    0    0    0    0    0    0
## 9   145    2   14   23    2    3    7    0   11
## 10    0    0    0    0    0    0    0    0    0
## 11    0    0    0    0    0    0    0    0    0
## 12    0    1    9    0    0    0    0    0    0
## 13    0    0    0    0    0    0    0    0    0
## 14    0    0    0    0    0    0    0    0    0
## 15    0    0    0    0    0    0    0    0    0
## 16    0    0    0    0    0    0    0    0    0
## 17    0    0    0    0    0    0    0    0    0

imm_gu <- geneUsage(immdata$data[c(1, 2)], "hs.trbv", .norm = T)
vis(imm_gu)
vis(imm_gu, .grid = T)
vis(imm_gu, .by = "Status", .meta = immdata$meta, .plot = "box")
imm_gu_cor <- geneUsageAnalysis(imm_gu, .method = "cor", .verbose = F)
vis(imm_gu_cor, .title = "Gene usage correlation", .leg.title = "Cor", .text.size = 1.5)
imm_cl_pca <- geneUsageAnalysis(imm_gu, "js+pca+kmeans", .verbose = F)
vis(imm_cl_pca, .plot = "clust")
imm_cl_pca2 <- geneUsageAnalysis(imm_gu, "js+pca+kmeans", .k = 3, .verbose = F)
vis(imm_cl_pca2)
p1 <- vis(spectratype(immdata$data[[1]], .quant = "id", .col = "nt"))
p2 <- vis(spectratype(immdata$data[[1]], .quant = "count", .col = "aa+v"))
p1 + p2

现在是不是觉得immunarch包很nice,同时值得提醒的是immunarch提供的可视化是基于ggplot2的vis函数来展示结果,这样对可视化做起来修改那也是相当的方便。该包还有更多的分析功能,我这里就不一一展示了,留个大家自己去探索吧。

最后

TCR-seq关键要理解这些数据结构以及其背后能够反应的生物学意义,这一点我也在继续学习中。。。




原文地址:访问原文地址
快照地址: 访问文章快照