setwd(‘D:/tcgabio‘)
rm(list = ls())

# TCGA-12-4567-01-blah-blah --> 這是Normal
# TCGA-12-4567-11-blah-blah --> 這是tumor
# 註意黑體的部分。01-09是tumor；10-19是Normal；20-29是Control
library(TCGAbiolinks)
# 下載前的query
query <- GDCquery(project = "TCGA-COAD", 
                    data.category = "Transcriptome Profiling", 
                    data.type = "Gene Expression Quantification", 
                    workflow.type = "HTSeq - FPKM-UQ")
GDCdownload(query)
# 將下載好的query轉換成一個SummerizedExperiment的文件，這個以rda為後綴的文件是一個總結性文件，
# 有了它，我們可以不再需要之前下載的raw數據，所以後面的remove.files.prepared可以選擇True，
# 這樣會把之前下載的大量文件刪除，當然也可以留著不刪除（即default）。
dataCOAD <- GDCprepare(query, save = TRUE, 
                       save.filename = "dataCOAD_summerizedExperiment.rda",
                       remove.files.prepared = TRUE)

# 可以看一看rda文件，用到的package是SummarizedExperiment
library(SummarizedExperiment)
samples.information=colData(dataCOAD)

# 數據準備好了，我們接下來開始進行DEA分析。所謂DEA，也就是Differential Expression Analysis，將Tumor組和對照組進行比較。
# 首先，將剛才GDCprepare好的數據進行normalization，用normalization()
# 這裏註意geneInfo=geneInfoHT，default其實是geneInfo，但由於我們前面選擇的是HTseq，所以要選擇geneInfoHT


dataNorm <- TCGAanalyze_Normalization(tabDF = dataCOAD, geneInfo = geneInfoHT)

# 之後，常規選擇，用Filtering()

dataFilt <- TCGAanalyze_Filtering(tabDF = dataNorm,
                                  method ="quantile", 
                                  qnt.cut = 0.25)

# 接著，定義對照組（這裏的對照組是Solid normal tissue），用到SampleType()，定義腫瘤組，用SampleType()

samplesNT <- TCGAquery_SampleTypes(barcode = colnames(dataFilt),
                                   typesample = c("NT"))

samplesTP <- TCGAquery_SampleTypes(barcode = colnames(dataFilt), 
                                   typesample = c("TP"))

# 進行DEA分析，用到DEA()

dataDEGs <- TCGAanalyze_DEA(mat1 =dataFilt[,samplesNT],
                            mat2 = dataFilt[,samplesTP],
                            Cond1type = "Normal",
                            Cond2type = "Tumor",
                            fdr.cut = 0.01 ,
                            logFC.cut = 1,
                            method = "glmLRT")
# 最後，將分析好的數據整入進一個表格裏，用到LevelTab()

dataDEGsFiltLevel <- TCGAanalyze_LevelTab(dataDEGs,"Tumor","Normal",
                                          dataFilt[,samplesTP],dataFilt[,samplesNT])

# 將表格保存到一個csv的文件
write.csv(dataDEGsFiltLevel,file="DEA_COAD.csv")

#最後得到得csv文件如下：

TCGAbiolinks（知乎整理）