immunarch — Fast and Seamless Exploration of Single-cell and Bulk T-cell/Antibody Immune Repertoires in R
10× Genomics單細胞免疫組庫VDJ分析必知必會
免疫組庫數(shù)據(jù)分析||immunarch教程：克隆型分析
免疫組庫數(shù)據(jù)分析||immunarch教程：探索性數(shù)據(jù)分析
免疫組庫數(shù)據(jù)分析||immunarch教程：載入10X數(shù)據(jù)
免疫組庫數(shù)據(jù)分析||immunarch教程：快速開始
免疫組庫數(shù)據(jù)分析||immunarch教程：GeneUsage分析
免疫組庫數(shù)據(jù)分析||immunarch教程：Diversity 分析
今天，我們繼續(xù)我們的免疫組庫數(shù)據(jù)分析的Demos，這一次我們來談?wù)凜lonotype tracking分析。像我這樣剛?cè)腴T免疫組庫的人首先會問什么是Clonotype tracking？
克隆型追蹤(Clonotype tracking)是監(jiān)測疫苗接種和癌癥免疫感興趣的克隆型頻率變化的常用方法。例如，研究人員可以在疫苗接種前和接種后的不同時間點跟蹤克隆類型，或分析腫瘤樣本中惡性克隆型的生長。所謂的追蹤用到的可視化方法類似我們提到過的?；鶊D，具體地：?；鶊D在單細胞數(shù)據(jù)探索中的應(yīng)用。immunarch中集成了多種克隆型跟蹤方法。目前有三種方法可供選擇。trackClonotypes的輸出可以立即用vis函數(shù)可視化。
追蹤最豐富的克隆型
最簡單的方法是從一個輸入免疫序列中選擇最豐富的克隆類型，并批量跟蹤所有的免疫序列。參數(shù).which和.col用于選擇免疫序列、從中獲取的克隆類型數(shù)量以及使用的列。從第一個庫中選擇5個最豐富的克隆型，并利用其CDR3核苷酸序列對其進行跟蹤:
library(immunarch); data(immdata) # Load the package and the test dataset
?trackClonotypes
tc1 <- trackClonotypes(immdata$data, list(1, 5), .col = "nt")
tc1
CDR3.nt
1: TGCGCCAGCAGCCAAGAAGGGACAGGGTATTCCGGGGAGCTGTTTTTT
2: TGCGCCAGCAGCTACAGGGTTGGCACAGATACGCAGTATTTT
3: TGTGCCACCAGCACCAACAGGGGCGGAACCCCAGCAGATACGCAGTATTTT
4: TGTGCCACCAGCATCGGAGGCGGGAGCTACGAGCAGTACTTC
5: TGTGCCAGCAGTCCTTGGACAGGGAGTATGGCCCTCCACTTT
A2-i129 A2-i131 A2-i133 A2-i132 A4-i191 A4-i192 MS1
1: 0.020352941 0 0 0 0 0 0
2: 0.019176471 0 0 0 0 0 0
3: 0.007764706 0 0 0 0 0 0
4: 0.006352941 0 0 0 0 0 0
5: 0.005647059 0 0 0 0 0 0
MS2 MS3 MS4 MS5 MS6
1: 0 0 0 0 0
2: 0 0 0 0 0
3: 0 0 0 0 0
4: 0 0 0 0 0
5: 0 0 0 0 0
which參數(shù)(第二個參數(shù))的值list(1,5)意味著從輸入字符串immdata$data的第一個列表中選擇5個clonotypes。col參數(shù)的值“nt”意味著函數(shù)應(yīng)該只接受CDR3核苷酸序列。
從“MS1”庫中選擇10個最豐富的氨基酸克隆型序列及其V基因進行跟蹤:
tc2 <- trackClonotypes(immdata$data, list("MS1", 10), .col = "aa+v")
tc2
CDR3.aa V.name A2-i129 A2-i131 A2-i133
1: CASSFEGAMDTQYF TRBV7-6 0 0 0
2: CASSLGDSTYEQYF TRBV5-6 0 0 0
3: CASSLGLREQGETQYF TRBV28 0 0 0
4: CASSLQAGGNTDTQYF TRBV7-2 0 0 0
5: CASSLYSNEQFF TRBV7-9 0 0 0
6: CASSVYSTISEQYF TRBV9 0 0 0
7: CSARDLANSYEQYF TRBV20-1 0 0 0
8: CSTEEDSYNEQFF TRBV20-1 0 0 0
9: CSVELRTESGYEQYF TRBV29-1 0 0 0
10: CSYRTGGPEQYF TRBV29-1 0 0 0
A2-i132 A4-i191 A4-i192 MS1 MS2
1: 0 0 0 0.008941176 0.0000000000
2: 0 0 0 0.011529412 0.0000000000
3: 0 0 0 0.007058824 0.0001176471
4: 0 0 0 0.063529412 0.0000000000
5: 0 0 0 0.004470588 0.0000000000
6: 0 0 0 0.037647059 0.0000000000
7: 0 0 0 0.009529412 0.0000000000
8: 0 0 0 0.024000000 0.0000000000
9: 0 0 0 0.017764706 0.0000000000
10: 0 0 0 0.007294118 0.0000000000
MS3 MS4 MS5 MS6
1: 0.0000000000 0.0000000000 0 0
2: 0.0000000000 0.0001176471 0 0
3: 0.0000000000 0.0000000000 0 0
4: 0.0000000000 0.0000000000 0 0
5: 0.0000000000 0.0000000000 0 0
6: 0.0001176471 0.0001176471 0 0
7: 0.0000000000 0.0000000000 0 0
8: 0.0001176471 0.0000000000 0 0
9: 0.0000000000 0.0000000000 0 0
10: 0.0000000000 0.0000000000 0 0
參數(shù)的值list("MS1"， "10")，它的意思是選擇10個clonotypes從命名為"MS1"的指令集在指令集immdata$data的輸入列表中。col的“aa+v”值意味著該功能應(yīng)該同時取CDR3氨基酸序列和最豐富的克隆型的v基因片段。
同樣經(jīng)典而又傻瓜式地操作
p1 <- vis(tc1)
p2 <- vis(tc2)
p1 / p2
用特定的核苷酸或氨基酸序列追蹤克隆型
為了跟蹤特定的clonotype序列，可以提供核苷酸或氨基酸序列作為which參數(shù)，同時提供列.col，以指定在哪些列中搜索序列。例如，要跟蹤下面指定的七個CDR3氨基酸序列，你需要執(zhí)行以下代碼:
target <- c("CASSLEETQYF", "CASSDSSGGANEQFF", "CASSDSSGSTDTQYF", "CASSLAGGYNEQFF", "CASSDSAGGTDTQYF", "CASSLDSYEQYF", "CASSSAGGYNEQFF")
tc <- trackClonotypes(immdata$data, target, .col = "aa")
vis(tc)
利用特定序列和基因片段追蹤克隆型
與之前的方法相比，利用序列和基因片段的信息來追蹤克隆型是可能的。你有一個特定的CDR3序列和基因片段的數(shù)據(jù)框。我們將通過從批次的第一個清單中選擇10個最豐富的克隆類型來模擬這一點:
target <- immdata$data[[1]] %>%
select(CDR3.aa, V.name) %>%
head(10)
target
# A tibble: 10 x 2
CDR3.aa V.name
<chr> <chr>
1 CASSQEGTGYSGELFF TRBV4-1
2 CASSYRVGTDTQYF TRBV4-1
3 CATSTNRGGTPADTQYF TRBV15
4 CATSIGGGSYEQYF TRBV15
5 CASSPWTGSMALHF TRBV27
6 CASQGDSFNSPLHF TRBV4-1
7 CASSQDMGGRNTGELFF TRBV4-1
8 CASSEEPRLFGYTF TRBV2
9 CASSQPGQGGGDEQFF TRBV4-1
10 CASSWVARGPYEQYF TRBV6-6
tc <- trackClonotypes(immdata$data, target)
vis(tc)
請注意，您可以使用目標數(shù)據(jù)框中的任何列，例如CDR3核苷酸和氨基酸序列以及任何基因片段。所以如何確定哪一列呢？
可視化跟蹤
根據(jù)你的研究和審美需求，有三種可視化克隆型追蹤的方法。要選擇圖形的類型，需要提供“。- .plot = "smooth" -默認使用，使用光滑線條和堆疊條形圖進行可視化;- .plot = "area" -使用豐度線下面的區(qū)域來可視化豐度;- .plot =“l(fā)ine”-只可視化線，連接同一克隆型在時間點之間的豐度水平。
target <- c("CASSLEETQYF", "CASSDSSGGANEQFF", "CASSDSSGSTDTQYF", "CASSLAGGYNEQFF", "CASSDSAGGTDTQYF", "CASSLDSYEQYF", "CASSSAGGYNEQFF")
tc <- trackClonotypes(immdata$data, target, .col = "aa")
p1 <- vis(tc, .plot = "smooth")
p2 <- vis(tc, .plot = "area")
p3<- vis(tc, .plot = "line")
library(patchwork)
p1+ p2 +p3
改變樣本的順序
vis函數(shù)的order參數(shù)控制可視化中樣本的順序。您可以通過您計劃可視化的樣本索引或樣本名稱。
# Passing indices
names(immdata$data)[c(1, 3, 5)] # check sample names
vis(tc, .order = c(1, 3, 5))
# You can change the order
vis(tc, .order = c(5, 1, 3))
當然，這樣也是可以的，圖我就不貼了啊。
# Passing sample names
vis(tc, .order = c("A2-i129", "A2-i133", "A4-i191"))
如果元數(shù)據(jù)（metadata）包含有關(guān)時間的信息，如接種疫苗或腫瘤樣本的時間點，則可以使用它相應(yīng)地對樣本重新排序。在我們的示例中，immdata$meta不包含關(guān)于時間點的信息，因此我們將模擬這種情況。
immdata$meta$Timepoint <- sample(1:length(immdata$data))
immdata$meta
# A tibble: 12 x 7
Sample ID Sex Age Status Lane Timepoint
<chr> <chr> <chr> <dbl> <chr> <chr> <int>
1 A2-i129 C1 M 11 C A 2
2 A2-i131 C2 M 9 C A 11
3 A2-i133 C4 M 16 C A 7
4 A2-i132 C3 F 6 C A 3
5 A4-i191 C8 F 22 C B 1
6 A4-i192 C9 F 24 C B 4
7 MS1 MS1 M 12 MS C 10
8 MS2 MS2 M 30 MS C 6
9 MS3 MS3 M 8 MS C 12
10 MS4 MS4 F 14 MS C 8
11 MS5 MS5 F 15 MS C 5
12 MS6 MS6 F 15 MS C 9
sample_order <- order(immdata$meta$Timepoint)
immdata$meta$Timepoint[sample_order]
immdata$meta$Sample[sample_order]
vis(tc, .order = sample_order)
vis(tc, .order = order(immdata$meta$Timepoint))
改變調(diào)色板
在R控制臺中運行?scale_fill_brewer，了解ColorBrewer及其配色方案的更多信息
p1<- vis(tc) + scale_fill_brewer(palette = "Spectral")
p2 <- vis(tc) + scale_fill_brewer(palette = "RdBu")
p1 + p2
我們?yōu)槭裁匆芯緾lonotype tracking??？是為了看Clonotype 在不同樣本中的分布情況進一步看Clonotype的可能的遷移和轉(zhuǎn)化狀態(tài)，所以關(guān)鍵的是樣本分組的生物學意義啊。
