一、TCGA數據下載

1、TCGAbiolinks包

• 介紹：可用于檢索，下載，并準備TCGA數據用于下游分析的R包
• 優點：具備一體化的下載整合，無需再使用復雜的方法對下載的單個數據重新進行整合。
• 使用GDCquery_clinic(project , type = "clinical") //臨床數據下載 qurey <- GDCquery(project, data.category, data.type, workflow.type,legacy = FALSE, access, platform, file.type, barcode,experimental.strategy, sample.type) GDCdownload(query, method = "api", files.per.chunk = 100) data <- GDCprepare(query) 1. projectTCGAbiolinks:::getGDCprojects()$project_id //查看有哪些癌癥類型2. data.categoryTCGAbiolinks:::getProjectSummary(project) //查看某癌癥有哪些數據類型3. data.typedata.type = "Gene Expression Quantification"//下載rna-seq的counts數據data.type = "miRNA Expression Quantification"//下載miRNA數據data.type = "Copy Number Segment" //下載Copy Number Variation數據......4. workflow.typeworkflow.type="HTSeq - Counts" //原始count數workflow.type="HTSeq - FPKM-UQ" //FPKM上四分位數標準化值workflow.type="HTSeq - FPKM" //FPKM值/表達量值......5. legacy：不同的數據來源 Legacy 與 harmonizedGDC Legacy Archive:以前在CGHUB和TCGA數據門戶中存儲的數據的原始數據，由TCGA數據協調中心（DCC）托管，在該門戶中用GRCH37（HG19）和GRCH36（HG18）作為參考基因組 GDC harmonized database:可用數據與grch38（hg38）使用gdc生物信息學流程進行協調，該流程提供了生物標本和臨床數據標準化的方法，簡單講就是對數據進行了一定標準化處理。harmonized數據庫包括轉錄譜數據，甲基化數據，miRNA數據，但缺少芯片數據6. 其他參數根據需求選擇，可通過TCGA官網查看詳細內容

二、免疫浸潤

1、immunedeconv包

輸入要求：行名為基因名(HGNC symbols )，列名為樣本名的矩陣
使用：

res = deconvolute(gene_expression_matrix, method)

特殊情況

• CIBERSORT：//需要注明CIBERSORT.R與LM22.txt文件路徑set_cibersort_binary("/path/to/CIBERSORT.R")set_cibersort_mat("/path/to/LM22.txt")
• TIMER:deconvolute(your_mixture_matrix, "timer",indications=c("SKCM", "SKCM", "BLCA")) //indications是一個表示每個樣本癌癥類型的向量，TIMER支持的癌癥類型查看可 immunedeconv::timer_available_cancers
• xCell
  在cell_type_map中xCell 僅有39種細胞類型，故使用xCell網站運行：https://xcell.ucsf.edu/

三、生存分析

區別

• 統計描述
  估計生存率（生存函數）：Kaplan-Meier法

• 統計推斷
  比較各組的生存率（生存曲線）：log-rank檢驗——單因素/分層分析
  生存期長短的影響因素分析： COX回歸——單因素/多因素分析

K-M曲線可以很直觀地表現出兩組或多組的生存率或死亡率，因為生存數據往往都是非正態分布，因此常使用非參數的檢驗方法log-rank檢驗。在實際寫作中，K-M曲線都會搭配log-rank檢驗的P值。

Cox回歸本質上是一種回歸模型，它沒有直接使用生存時間，而是使用了危險度（hazard）作為因變量，該模型不用于估計生存率，而是用于因素分析，也就是找到某一個危險因素對結局事件發生的貢獻度。

因此可使用K-M方法繪制生存曲線，使用log-rank檢驗評估生存差異。使用Cox比例風險模型進行多因素分析，評估預后因素作用。

創建生存對象：Surv()
詳見 https://www.rdocumentation.org/packages/survival/versions/3.2-3/topics/Surv

//useage Surv(time, time2, event,type=c('right', 'left', 'interval', 'counting', 'interval2', 'mstate'),origin=0) //example Surv(OSTime,OS)

擬合生存曲線：survfit()
詳見https://www.rdocumentation.org/packages/survival/versions/3.2-3/topics/survfit

//usegae survfit(formula, ...) //example survfit(Surv(OSTime, OS) ~ variable, data = data)

繪制生存曲線：ggsurvplot()
詳見https://www.rdocumentation.org/packages/survminer/versions/0.4.7/topics/ggsurvplot

//useage ggsurvplot(fit, //survfit()對象/對象列表/包含生存曲線摘要的數據框data = NULL, //用于擬合生存曲線的數據集，如未提供則從fit中提取fun = NULL,color = NULL, //默認情況下，使用參數color =“ strata”按層為生存曲線著色，也可以使用參數調色板指定自定義調色板palette = NULL, //調色板linetype = 1, //線型conf.int = FALSE, //繪制置信區間pval = FALSE, //繪制P值pval.method = FALSE, //Pval為T時，顯示P值檢驗方法test.for.trend = FALSE,surv.median.line = "none",risk.table = FALSE,cumevents = FALSE,cumcensor = FALSE,tables.height = 0.25,group.by = NULL,facet.by = NULL,add.all = FALSE,combine = FALSE,ggtheme = theme_survminer(),tables.theme = ggtheme,... )

多因素分析：coxph()
詳見：https://www.rdocumentation.org/packages/survival/versions/3.2-3/topics/coxph

coxph(formula, //~左側surv函數返回的生存對象，右側為待分析因素data=, //用于解釋formula/weights/subset中命名變量的數據框weights, subset, na.action, init, control, ties=c("efron","breslow","exact"), //默認為efronsingular.ok=TRUE, robust, model=FALSE, x=FALSE, y=TRUE, tt, method=ties,id, cluster, istate, statedata, ...)

四、 差異分析

1、使用FPKM數據進行差異基因分析

//將FPKM數據轉換為TPM數據 expMatrix <- FPKM_matrix fpkmToTpm <- function(fpkm) {exp(log(fpkm) - log(sum(fpkm)) + log(1e6)) }apply(expMatrix,2,fpkmToTpm)//限定分類順序 Group <- factor(group,levels = ,ordered = F) design2 <- model.matrix(~factor( Group ))//表達矩陣數據校正 exprSet <-data boxplot(exprSet,outline=FALSE, notch=T,col=Group, las=2)//查看數據是否需要標準化 library(limma) exprSet=normalizeBetweenArrays(exprSet)//數據標準化 boxplot(exprSet,outline=FALSE, notch=T,col=Group, las=2) exprSet <- log2(exprSet+1)//需判斷數據是否要log化//差異分析： dat <- exprSet fit=lmFit(dat,design2) fit=eBayes(fit) options(digits = 4) deg=topTable(fit,coef=2,adjust='BH',number = Inf) difflimma <- subset(deg,abs(deg$logFC)>1&deg$P.Value<0.05) //篩選條件為|logFC|>1，Pvalue<0.05

五、富集分析

library(clusterProfiler) library(org.Hs.eg.db) //id轉換 gene_id =bitr(gene_id,fromType = "SYMBOL",toType = "ENTREZID", OrgDb="org.Hs.eg.db")

1、GO分析

詳見https://www.rdocumentation.org/packages/clusterProfiler/versions/3.0.4/topics/enrichGO

enrichGO(gene, //entrez idOrgDb, //組織數據庫keytype = "ENTREZID", ont = "MF", //MF", "BP", "CC"，"ALL"pvalueCutoff = 0.05, pAdjustMethod = "BH", universe, qvalueCutoff = 0.2, minGSSize = 10, maxGSSize = 500, readable = FALSE)

2、KEGG分析

詳見https://www.rdocumentation.org/packages/clusterProfiler/versions/3.0.4/topics/enrichKEGG

enrichKEGG(gene, organism = "hsa", keyType = "kegg", //"kegg", 'ncbi-geneid', 'ncib-proteinid','uniprot'pvalueCutoff = 0.05, pAdjustMethod = "BH", universe, minGSSize = 10, maxGSSize = 500, qvalueCutoff = 0.2, use_internal_data = FALSE)

六、ggplot

元素

plot：整張圖，包括background和title

axis：包括stick，text，title

legend：包括background，text，title

facet：可分為外部stript（包括backgroud和text）和內部panel部分（包括backgroud、boder和網格線grid，其中粗的叫grid.major，細的叫grid.minor）。

函數類型

用于運算（如fortify_，mean_等）

初始化、展示（ggplot，plot，print等）

按變量組圖（facet_等）

繪圖（stat_，geom_，annotate），核心函數。

微調圖型：
- scale_：與aes內的各種美學（shape、color、fill、alpha）調整有關的函數。
- guides：調整所有的text。
- coord_：調整坐標。
- theme：調整不與數據有關的圖的元素的函數。

1、ggpubr

介紹：基于ggplot2的用于繪制符合出版物要求的圖形的可視化包
鏈接：https://www.jianshu.com/p/678213d605a5

分布圖(Distribution)

#構建數據集 set.seed(1234) df <- data.frame( sex=factor(rep(c("f", "M"), each=200)), weight=c(rnorm(200, 55), rnorm(200, 58))) head(df) ## sex weight ## 1 f 53.79293 ## 2 f 55.27743 ## 3 f 56.08444 ## 4 f 52.65430 ## 5 f 55.42912 ## 6 f 55.50606

密度分布圖以及邊際地毯線并添加平均值線

ggdensity(df, x="weight", add = "mean", rug = TRUE, color = "sex", fill = "sex", palette = c("#00AFBB", "#E7B800"))

帶有均值線和邊際地毯線的直方圖

gghistogram(df, x="weight", add = "mean", rug = TRUE, color = "sex", fill = "sex", palette = c("#00AFBB", "#E7B800"))

箱線圖與小提琴圖

#加載數據集ToothGrowth data("ToothGrowth") df1 <- ToothGrowth head(df1) ## len supp dose ## 1 4.2 VC 0.5 ## 2 11.5 VC 0.5 ## 3 7.3 VC 0.5 ## 4 5.8 VC 0.5 ## 5 6.4 VC 0.5 ## 6 10.0 VC 0.5 p <- ggboxplot(df1, x="dose", y="len", color = "dose", palette = c("#00AFBB", "#E7B800", "#FC4E07"), add = "jitter", shape="dose")#增加了jitter點，點shape由dose映射p

增加不同組間的p-value值，可以自定義需要標注的組間比較

my_comparisons <- list(c("0.5", "1"), c("1", "2"), c("0.5", "2")) p+stat_compare_means(comparisons = my_comparisons)+#不同組間的比較 stat_compare_means(label.y = 50)

內有箱線圖的小提琴圖

ggviolin(df1, x="dose", y="len", fill = "dose", palette = c("#00AFBB", "#E7B800", "#FC4E07"), add = "boxplot", add.params = list(fill="white"))+ stat_compare_means(comparisons = my_comparisons, label = "p.signif")+#label這里表示選擇顯著性標記（星號） stat_compare_means(label.y = 50)

條形圖

data("mtcars") df2 <- mtcars df2$cyl <- factor(df2$cyl) df2$name <- rownames(df2)#添加一行name head(df2[, c("name", "wt", "mpg", "cyl")])

按從小到大順序繪制條形圖

不分組排序 ggbarplot(df2, x="name", y="mpg", fill = "cyl", color = "white", palette = "jco",#雜志jco的配色 sort.val = "desc",#下降排序 sort.by.groups=FALSE,#不按組排序 x.text.angle=60) 按組進行排序 ggbarplot(df2, x="name", y="mpg", fill = "cyl", color = "white", palette = "jco",#雜志jco的配色 sort.val = "asc",#上升排序,區別于desc，具體看圖演示 sort.by.groups=TRUE,#按組排序 x.text.angle=90)

偏差圖
偏差圖展示了與參考值之間的偏差

df2$mpg_z <- (df2$mpg-mean(df2$mpg))/sd(df2$mpg) df2$mpg_grp <- factor(ifelse(df2$mpg_z<0, "low", "high"), levels = c("low", "high")) head(df2[, c("name", "wt", "mpg", "mpg_grp", "cyl")])

繪制排序過的條形圖

ggbarplot(df2, x="name", y="mpg_z", fill = "mpg_grp", color = "white", palette = "jco", sort.val = "asc", sort.by.groups = FALSE, x.text.angle=60, ylab = "MPG z-score", xlab = FALSE, legend.title="MPG Group")

坐標軸變換

ggbarplot(df2, x="name", y="mpg_z", fill = "mpg_grp", color = "white", palette = "jco", sort.val = "desc", sort.by.groups = FALSE, x.text.angle=90, ylab = "MPG z-score", xlab = FALSE, legend.title="MPG Group", rotate=TRUE, ggtheme = theme_minimal())

點圖(Dot charts)
棒棒糖圖(Lollipop chart)
棒棒圖可以代替條形圖展示數據

ggdotchart(df2, x="name", y="mpg", color = "cyl", palette = c("#00AFBB", "#E7B800", "#FC4E07"), sorting = "ascending", add = "segments", ggtheme = theme_pubr())

可以自設置各種參數

ggdotchart(df2, x="name", y="mpg", color = "cyl", palette = c("#00AFBB", "#E7B800", "#FC4E07"), sorting = "descending", add = "segments", rotate = TRUE, group = "cyl", dot.size = 6, label = round(df2$mpg), font.label = list(color="white", size=9, vjust=0.5), ggtheme = theme_pubr())

偏差圖

ggdotchart(df2, x="name", y="mpg_z", color = "cyl", palette = c("#00AFBB", "#E7B800", "#FC4E07"), sorting = "descending", add = "segment", add.params = list(color="lightgray", size=2), group = "cyl", dot.size = 6, label = round(df2$mpg_z, 1), font.label = list(color="white", size=9, vjust=0.5), ggtheme = theme_pubr())+ geom_line(yintercept=0, linetype=2, color="lightgray")

Cleveland點圖

ggdotchart(df2, x="name", y="mpg", color = "cyl", palette = c("#00AFBB", "#E7B800", "#FC4E07"), sorting = "descending", rotate = TRUE, dot.size = 2, y.text.col=TRUE, ggtheme = theme_pubr())+ theme_cleveland()

2、色板擴充

默認顏色brewer.pal.info #查詢默認顏色

• seq類型：單漸變色，一種主色由淺到深
• qual類型：區分色，幾種區分度很高的顏色組合
• div類型：雙漸變色，一種顏色到另外一種顏色的漸變，有兩種主色

色板擴充//直接使用默認色板 ggplot(...)+scale_fill_brewer(palette="Set1")+...//colorRamPalette使用默認色板 library(RColorBrewer) colors <-colorRampPalette(brewer.pal(9,"Set1")) //colorRamPalette使用自定義色板 colors <-colorRampPalette(c("BlueViolet","MediumOrchid","Orchid","MediumSlateBlue","CornflowerBlue","LightSteelBlue","LightSkyBlue","SkyBlue","PowDerBlue","PaleTurquoise","LightSeaGreen","Turquoise","LightGreen","DarkSeaGreen","Khaki","NavajoWhite","Gold","Orange")) ggplot(...+...+ scale_fill_manual(values = colors(22))+...)//22表示漸變梯度

總結

以上是生活随笔為你收集整理的2020.07.01-07.15学习小结的全部內容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網站內容還不錯，歡迎將生活随笔推薦給好友。