引言：在前面我們了解了如何使用TCGAbiolinks檢索并獲取TCGA數(shù)據(jù)庫的公開數(shù)據(jù)。今天小編就用前面涉及到的代碼，下載今天數(shù)據(jù)準備需要用到的TCGA樣本數(shù)據(jù)。

一、數(shù)據(jù)下載階段

第一步：GDCquery（）篩選我們需要的數(shù)據(jù)，TCGAbiolinks包下載TCGA數(shù)據(jù)進行表達差異分析－肝癌案例

library（＂TCGAbiolinks＂）

query ＜－ GDCquery（project ＝ ＂TCGA－LIHC＂，

data．category ＝ ＂Transcriptome Profiling＂，

data．type ＝ ＂Gene Expression Quantification＂，

workflow．type ＝ ＂HTSeq － Counts＂）

上圖為通過TCGA GDC鏈接中根據(jù)篩選條件查看的符合要求結(jié)果。下圖為通過GDCquery（）函數(shù)中傳入對應(yīng)的參數(shù)得到的結(jié)果。兩者對比，我們可以發(fā)現(xiàn)，兩者是一模一樣的。說明代碼執(zhí)行正確。前面一期中，我們有詳細談及 GDCquery，可做參考。

samplesDown ＜－ getResults（query，cols＝c（＂cases＂））

＃getResults（query， rows， cols）根據(jù)指定行名或列名從query中獲取結(jié)果，此處用來獲得樣本的barcode

＃ 此處共檢索出424個barcodes

getResults（）中用到的參數(shù)：

參數(shù)用法query

來自GDCquery的結(jié)果rows用于指定特定的行cols用于指定特定的列

＃ 從samplesDown中篩選出TP（實體腫瘤）樣本的barcodes

＃ TCGAquery＿SampleTypes（barcode， typesample）

＃ TP代表PRIMARY SOLID TUMOR；NT－代表Solid Tissue Normal（其他組織樣本可參考學習文檔）

＃＃此處共檢索出371個TP樣本barcodes

dataSmTP ＜－ TCGAquery＿SampleTypes（barcode ＝ samplesDown，

typesample ＝ ＂TP＂）

＃ 從samplesDown中篩選出NT（正常組織）樣本的barcode

＃此處共檢索出50個NT樣本barcodes

dataSmNT ＜－ TCGAquery＿SampleTypes（barcode ＝ samplesDown，

typesample ＝ ＂NT＂）

TCGAquery＿SampleTypes中的參數(shù)詳解：

參數(shù)用法barcodeTCGA中的barcodes列表typesample用于指定篩選哪種類型的組織樣本，如腫瘤組織“TP”，正常組織“NT”

補充TCGA中的組織樣本類型：

TPPRIMARY SOLID TUMORTMMetastaticTRRECURRENT SOLID TUMORTAMAdditional MetastaticTBPrimary Blood Derived Cancer－Peripheral BloodTHOCHuman Tumor Original CellsTRBMRecurrent Blood Derived Cancer－Bone MarrowTBM Primary Blood Derived Cancer－Bone MarrowTAPAdditional－New PrimaryNB Blood Derived Normal NTSolid Tissue NormalNBCBuccal Cell Normal？？？NEBVEBV Immortalized NormalNBMBone Marrow Normal

＃＃＃設(shè)置barcodes參數(shù)，篩選符合要求的371個腫瘤樣本數(shù)據(jù)和50正常組織數(shù)據(jù)

queryDown ＜－ GDCquery（project ＝ ＂TCGA－LIHC＂，

data．category ＝ ＂Transcriptome Profiling＂，

data．type ＝ ＂Gene Expression Quantification＂，

workflow．type ＝ ＂HTSeq － Counts＂，

barcode ＝ c（dataSmTP， dataSmNT））

＃barcode參數(shù)：根據(jù)傳入barcodes進行數(shù)據(jù)過濾

上圖為 queryDown＜－GDCquery（）的結(jié)果，僅選擇了選擇371個正常組織和50個腫瘤組織樣本。

第二步：GDCdownload（）下載GDCquery（）得到的結(jié)果

＃ 下載數(shù)據(jù)，默認存放位置為當前工作目錄下的GDCdata文件夾中。

GDCdownload（queryDown，method ＝ ＂api＂， directory ＝ ＂GDCdata＂，

files．per．chunk ＝ 10）

＃method ；＂API＂或者＂client＂。＂API＂速度更快，但是容易下載中斷。

＃directory：下載文件的保存地址。Default： GDCdata。

＃files．per．chunk ＝ NULL：使用API下載大文件的時候，可以把文件分成幾個小文件來下載，可以解決下載容易中斷的問題。

GDCdownload（query ＝ queryDown）

說明：由于小編前面已經(jīng)下載過該TCGA數(shù)據(jù)，所以這里顯示的是421個文件已存在。如果還沒有下載的話，可能需要根據(jù)自己的網(wǎng)速等待一些時間。

顯示這樣的結(jié)果，就算下載成功啦！文件默認保存在 Rstudio默認路徑下的GDCdata中。前面就是我們利用第一期知識進行數(shù)據(jù)下載環(huán)節(jié)，權(quán)當溫習功課吧——接下來我們就開始此期的數(shù)據(jù)處理～～

二、數(shù)據(jù)處理

第三步：GDCprepare（）將前面GDCquery（）的結(jié)果準備成R語言可處理的SE（SummarizedExperiment）文件。

＃讀取下載的數(shù)據(jù)并將其準備到R對象中，在工作目錄生成（save＝TRUE）LIHC＿case．rda文件

＃ GDCprepare（）：Prepare GDC data，準備GDC數(shù)據(jù)，使其可用于R語言中進行分析

dataPrep1 ＜－ GDCprepare（query ＝ queryDown， save ＝ TRUE， save．filename ＝

＂LIHC＿case．rda＂）

GDCprepare（）中的參數(shù)：

參數(shù)用法query來自GDCquery的結(jié)果save是否將結(jié)果保存為RData object，默認為TRUEsave．filename文件名，如果沒有設(shè)置，系統(tǒng)將默認設(shè)置directory文件數(shù)據(jù)的文件夾，默認為“GDCdata”summarizedExperiment是否生成summarizedExperiment對象，默認TRUE

第四步：TCGAanalyze＿Preprocessing（）對數(shù)據(jù)進行預處理：使用spearman相關(guān)系數(shù)去除數(shù)據(jù)中的異常值

＃ 去除dataPrep1中的異常值，dataPrep1數(shù)據(jù)中含有腫瘤組織和正常組織的數(shù)據(jù)

＃ TCGAanalyze＿Preprocessing（object， cor．cut ＝ 0， filename ＝ NULL，

width ＝ 1000， height ＝ 1000， datatype ＝ names（assays（object））［1］）

＃ 函數(shù)功能描述：Array Array Intensity correlation （AAIC） and correlation boxplot to define outlier

dataPrep2 ＜－ TCGAanalyze＿Preprocessing（object ＝ dataPrep1，

cor．cut ＝ 0．6，

datatype ＝ ＂HTSeq － Counts＂）

＃將預處理后的數(shù)據(jù)dataPrep2，寫入新文件“LIHC＿dataPrep．csv”

write．csv（dataPrep2，file ＝ ＂LIHC＿dataPrep．csv＂，quote ＝ FALSE）

這里將生成一個array－array intensity correlation（AAIC）相關(guān)性熱圖，如下：

TCGAanalyze＿Preprocessing（）中的參數(shù)：

參數(shù)用法object來自TCGAprepare的結(jié)果cor．cut設(shè)置閾值，根據(jù)樣本中各個樣本之間的spearman相關(guān)系數(shù)進行過濾。默認為0filename設(shè)置生成圖片文件的名稱，默認為PreprocessingOutput．pngwidth生成圖片的寬度？？ height生成圖片的高度datatype描述RangedSummarizedExperiment 數(shù)據(jù)類型的字符串

第五步：TCGAtumor＿purity（）篩選腫瘤純度大于60％的腫瘤barcodes

＃ TCGAtumor＿purity（barcodes， estimate， absolute， lump， ihc， cpe），使用來自5種方法的5個估計值作為閾值對TCGA樣本進行過濾，這5個值是estimate， absolute， lump， ihc， cpe，這里設(shè)置cpe＝0．6（cpe是派生的共識度量，是將所有方法的標準含量歸一化后的均值純度水平，以使它們具有相等的均值和標準差）

＃篩選腫瘤純度大于等于60％的樣本數(shù)據(jù)

purityDATA ＜－ TCGAtumor＿purity（colnames（dataPrep1）， 0， 0， 0， 0， 0．6）

＃ filtered 為被過濾的數(shù)據(jù)， pure＿barcodes是我們要的腫瘤數(shù)據(jù)

Purity．LIHC＜－purityDATA＄pure＿barcodes

normal．LIHC＜－purityDATA＄filtered

filtered 為被過濾的數(shù)據(jù)（為正常組織的數(shù)據(jù)barcodes）， pure＿barcodes是我們要的腫瘤樣本barcodes。

第六步：將腫瘤表達矩陣與正常組織表達矩陣合并，進行基因注釋

＃獲取腫瘤純度大于60％的340個腫瘤組織樣本＋50個正常組織樣本，共計390個樣本

puried＿data ＜－dataPrep2［，c（Purity．LIHC，normal．LIHC）］

第七步：進行表達矩陣基因注釋

?；蜃⑨專枰虞d“SummarizedExperiment”包，“SummarizedExperiment container”每個由數(shù)字或其他模式的類似矩陣的對象表示。行通常表示感興趣的基因組范圍和列代表樣品。

＃if （！requireNamespace（＂BiocManager＂， quietly ＝ TRUE））

install．packages（＂BiocManager＂）

＃BiocManager：：install（＂SummarizedExperiment＂） ＃沒有的需要執(zhí)行下載代碼

library（＂SummarizedExperiment＂）

rowData（dataPrep1） ＃傳入數(shù)據(jù)dataPrep1必須為SummarizedExperiment對象

＃ DataFrame with 56512 rows and 3 columns

＃ ensembl＿gene＿id external＿gene＿name original＿ensembl＿gene＿id

＃ ＜character＞ ＜character＞ ＜character＞

＃ ENSG00000000003 ENSG00000000003 TSPAN6 ENSG00000000003．13

＃ ENSG00000000005 ENSG00000000005 TNMD ENSG00000000005．5

＃ ENSG00000000419 ENSG00000000419 DPM1 ENSG00000000419．11

＃ ENSG00000000457 ENSG00000000457 SCYL3 ENSG00000000457．12

＃將結(jié)果寫入文件“puried．LIHC．cancer．csv”

rownames（puried＿data）＜－rowData（dataPrep1）＄external＿gene＿name

write．csv（puried＿data，file ＝ ＂puried．LIHC．csv＂，quote ＝ FALSE）

第八步：進行表達矩陣標準化和過濾，得到用于差異分析的表達矩陣

｀TCGAanalyze＿Normalization（）｀使用EDASeq軟件包標準化mRNA轉(zhuǎn)錄本和miRNA。

＃TCGAanalyze＿Normalization（）執(zhí)行EDASeq包中的如下功能：

1． EDASeq：：newSeqExpressionSet

2． EDASeq：：withinLaneNormalization

3． EDASeq：：betweenLaneNormalization

4． EDASeq：：counts

dataNorm ＜－ TCGAanalyze＿Normalization（tabDF ＝ puried＿data，

geneInfo ＝ geneInfo，

method ＝ ＂gcContent＂）

TCGAanalyze＿Normalization中的參數(shù)：

參數(shù)用法tabDFRNAseq表達矩陣，行代表基因，列代表樣本geneInfo關(guān)于geneLength和gcContent的20531個基因的矩陣，“geneInfoHT”和“geneInfo”可選。method選擇標準化的方法，基于’gcContent’ 或 ’geneLength’的標準化方法可選

＃將標準化后的數(shù)據(jù)再過濾，去除掉表達量較低（count較低）的基因，得到最終的數(shù)據(jù)

dataFilt ＜－ TCGAanalyze＿Filtering（tabDF ＝ dataNorm，

method ＝ ＂quantile＂，

qnt．cut ＝ 0．25）

str（dataFilt）

＃num ［1：13083， 1：340］ 274 2432 60347 1012 1947 ．．．

＃－ attr（＊， ＂dimnames＂）＝List of 2

＃ ．．＄ ： chr ［1：13083］ ＂A1BG＂ ＂A1CF＂ ＂A2M＂ ＂A4GALT＂ ．．．

＃ ．．＄ ： chr ［1：390］ ＂TCGA－DD－AAD5－01A－11R－A41C－07＂ ＂TCGA－DD－A4NO－01A－11R－A28V－07＂ ＂TCGA－EP－A2KA－01A－11R－A180－07＂ ＂TCGA－DD－AACP－01A－11R－A41C－07＂ ．．．

TCGAanalyze＿Filtering（）中的參數(shù)：

參數(shù)用法tabDF數(shù)據(jù)框或者矩陣，行代表基因，列代表來自TCGA的樣本method用于過濾較低count數(shù)的基因的方法，有’quantile’， ’varFilter’， ’filter1’， ’filter2’qnt．cut選擇均值作為過濾的閾值

最后將過濾后的數(shù)據(jù)寫入文件“TCGA＿LIHC＿final．csv”，就得到我們用于后續(xù)差異分析的表達文件：

write．csv（dataFilt，file ＝ ＂TCGA＿LIHC＿final．csv＂，quote ＝ FALSE）

＃保留的是390個樣本（前340腫瘤，后50正常組織）

今天的數(shù)據(jù)預處理就講到這里，接下來我們將分享：數(shù)據(jù)分析（差異表達分析、富集分析和聚類分析等）。如果你喜歡的話，就加入我們一起挖數(shù)據(jù)吧～～