R语言进行TCGA配对样本差异基因分析

已于 2024-04-22 16:23:15 修改
阅读量1.7w 收藏 156
点赞数 14
文章标签： python 机器学习 数据分析 人工智能 数据可视化
于 2021-03-30 22:19:00 首次发布
版权声明：本文为博主原创文章，遵循 CC 4.0 BY-SA 版权协议，转载请附上原文出处链接和本声明。
本文链接：https://blog.csdn.net/qazplm12_3/article/details/115339797
版权

之前的一个推文是从UCSC XENA获取TCGA的表达和表型数据，然后利用代码对表达数据进行了ID注释，以及mRNA、lncRNA和miRNA的区分筛选，最后将患者ID和临床信息进行配比，用于后续分析。详细内容见推文：《利用R代码从UCSC XENA下载mRNA, lncRNA, miRNA表达数据并匹配临床信息》。

本期我将继续上次的内容，从TCGA 546个头颈癌数据集（Tumor = 502，Normal = 44）中，提取出43对癌和癌旁样本，并做配对差异分析， 然后绘制某个基因的配对箱式图。实际上TCGA好多癌症比如头颈癌、肝癌等，都有癌和癌旁的数据，且癌和癌旁都是一一配对的关系。所以在分析癌vs.癌旁的过程中，可以选择普通的差异分析，例如头颈癌的502个癌vs.44癌旁；另外一种思路是从中挑选出配对的癌vs.癌旁进行配对差异分析，例如头颈癌的43个癌vs.43个癌旁。实际上，同一个数据的非配对分析和配对分析差异还是很大的，详见我之前写过的一个帖子《差异分析|DESeq2完成配对样本的差异分析》。

关于DEseq2配对差异分析，很少有帖子涉及 (生信宝典注：不是很少有帖子涉及，而是看有没有发现，配对的个体信息可以视作一个批次因素，在DESeq2差异基因分析和批次效应移除,   高通量数据中批次效应的鉴定和处理 - 系列总结和更新,    典型医学设计实验GEO数据分析 (step-by-step) - Limma差异分析、火山图、功能富集中都有类似的处理方式)。

总结来说，拿到配对设计的数据，如果不进行配对分析而用常规的差异分析，这样的结果可能会大不相同。因此，在分析数据的时候，一定要明白实验设计。

下面代码展示了如何从546个样本中挑选出一一配对的癌和癌旁数据，并进行DEseq2配对差异分析。

remove(list = ls())
##加载包
suppressMessages(library(DESeq2))
suppressMessages(library(dplyr))
library(ggplot2)
library(ggpubr)
library(ggthemes)

1. 加载数据

首先，读入上期推文中处理好的表达矩阵（仅含mRNA）和表型数据。

# 1.1 表达矩阵
expr_data <- read.csv("./Rawdata/TCGA_mRNA_count_log2.csv",header = T,row.names = 1,check.names=F)
expr_data[1:3,1:4]

##       TCGA.BB.4224.01 TCGA.H7.7774.01 TCGA.CV.6943.01 TCGA.CN.5374.01
## A4GNT        2.584963         0.00000         1.00000        2.321928
## AAAS        10.952741        10.43463        11.27671       10.750707
## AACS        11.098032        13.62822        12.52943       11.171802

# 1.2 表型数据
phenotype <- read.csv("./Rawdata/TCGA_HNSC_phenotype.csv",header = T,row.names = 1,check.names=F, stringsAsFactors = T)
phenotype$subject <- as.factor(substring(rownames(phenotype),1,12))
phenotype[1:3,c("group","subject")]

##                 group      subject
## TCGA.BB.4224.01 Tumor TCGA.BB.4224
## TCGA.H7.7774.01 Tumor TCGA.H7.7774
## TCGA.CV.6943.01 Tumor TCGA.CV.6943

# 检查一下表达数据和表型数据是否是一一对应
table(colnames(expr_data) == rownames(phenotype)) #T

## 
## TRUE 
##  546

表达数据和表型数据如果不是一一对应的话，运行下面这行代码:

# 利用match函数进行匹配
phenotype <- phenotype[match(colnames(expr_data),rownames(phenotype)),]
identical(colnames(expr_data),rownames(phenotype))

## [1] TRUE

2. 筛选出配对的样本

2.1 将546个样本分为两部分，Tumor组和Normal组

# Tumor组
Tumor_group <- phenotype[which(phenotype$group == "Tumor"),] #TCGA的Tumor样本,502个
table(Tumor_group$group) #502个癌

## 
## Nontumor    Tumor 
##        0      502

# Normal组
Normal_group <- phenotype[which(phenotype$group == "Nontumor"),] #TCGA的control样本,44
table(Normal_group$group) #44个癌旁

## 
## Nontumor    Tumor 
##       44        0

# 找到癌和癌旁能够进行一一配对的样本的信息
Normal_group <- Normal_group[Normal_group$subject %in% Tumor_group$subject,]
#仅保留配对的对照组
table(Normal_group$group) # 43个存在配对的癌旁

## 
## Nontumor    Tumor 
##       43        0

上面这行代码就获取了43个配对的癌旁样本信息，根据这个信息可以获得所有配对的癌+癌旁，代码如下：

# 获取表达矩阵
paired_expr_data <- expr_data[,substring(colnames(expr_data),1,12) %in% 
                                     substring(rownames(Normal_group),1,12)] 
# 实际上，这里的substring(rownames(Normal_group),1,12)和Normal_group$ID是一样的
paired_phe <- phenotype[match(colnames(paired_expr_data),row.names(phenotype)),]
paired_phe <- droplevels(paired_phe)
dim(paired_expr_data) # 获得了86个样本表达信息，其中43个癌，43个癌旁，是一一配对的关系

## [1] 18192    86

write.csv(paired_expr_data,"./Output/TCGA_HNSC_mRNA_count_paired_43vs43.csv")

然后，进行配对的DEseq2差异分析。这里只保留了至少在75%的样本中都有表达的基因。

3. DEseq2配对差异分析

3.1 配对的表达矩阵

# UCSC Xena下载的数据是经过log2（count+1）处理的
data = (2^paired_expr_data-1) %>% apply(2, as.integer)# 还原log2的处理
row.names(data) = row.names(paired_expr_data)

# 保留至少在75%的样本中都有表达的基因
keep_data <- rowSums(data> 0) >= floor(0.75*ncol(data))
table(keep_data)

## keep_data
## FALSE  TRUE 
##  2599 15593

data <- data[keep_data,]
data[1:4,1:4]

##       TCGA.CV.6943.01 TCGA.CV.6959.01 TCGA.CV.7438.11 TCGA.CV.7242.11
## AAAS             2479            3304            1900            2190
## AACS             5911            3881            4759            3070
## AADAC               7               5             449            2618
## AADAT             140              77             160             127

dim(data)

## [1] 15593    86

3.2 对应的表型数据和配对信息

table(paired_phe$subject,paired_phe$group)

##               
##                Nontumor Tumor
##   TCGA.CV.6933        1     1
##   TCGA.CV.6934        1     1
##   TCGA.CV.6935        1     1
##   TCGA.CV.6936        1     1
##   TCGA.CV.6938        1     1
##   TCGA.CV.6939        1     1
##   TCGA.CV.6943        1     1
##   TCGA.CV.6955        1     1
##   TCGA.CV.6956        1     1
##   TCGA.CV.6959        1     1
##   TCGA.CV.6960        1     1
##   TCGA.CV.6961        1     1
##   TCGA.CV.6962        1     1
##   TCGA.CV.7091        1     1
##   TCGA.CV.7097        1     1
##   TCGA.CV.7101        1     1
##   TCGA.CV.7103        1     1
##   TCGA.CV.7177        1     1
##   TCGA.CV.7178        1     1
##   TCGA.CV.7183        1     1
##   TCGA.CV.7235        1     1
##   TCGA.CV.7238        1     1
##   TCGA.CV.7242        1     1
##   TCGA.CV.7245        1     1
##   TCGA.CV.7250        1     1
##   TCGA.CV.7252        1     1
##   TCGA.CV.7255        1     1
##   TCGA.CV.7261        1     1
##   TCGA.CV.7406        1     1
##   TCGA.CV.7416        1     1
##   TCGA.CV.7423        1     1
##   TCGA.CV.7424        1     1
##   TCGA.CV.7425        1     1
##   TCGA.CV.7432        1     1
##   TCGA.CV.7434        1     1
##   TCGA.CV.7437        1     1
##   TCGA.CV.7438        1     1
##   TCGA.CV.7440        1     1
##   TCGA.H7.A6C4        1     1
##   TCGA.HD.8635        1     1
##   TCGA.HD.A6HZ        1     1
##   TCGA.HD.A6I0        1     1
##   TCGA.WA.A7GZ        1     1

3.3 配对差异分析

# DEseq2的配对差异分析跑的有点慢，43对我的电脑需要约10分钟
# 所以，这里我为了加快运行速度，选了5对5.
data = data [,order(paired_phe$subject)]
data = as.matrix(data[,1:10])
data[1:3,1:4]

##       TCGA.CV.6933.01 TCGA.CV.6933.11 TCGA.CV.6934.01 TCGA.CV.6934.11
## AAAS             2993            1328            2816            1977
## AACS             5124            4997            3522            4108
## AADAC              98               1              16              63

coldata_paired <- paired_phe[colnames(data),]
coldata_paired[,c('subject','group')]

##                      subject    group
## TCGA.CV.6933.01 TCGA.CV.6933    Tumor
## TCGA.CV.6933.11 TCGA.CV.6933 Nontumor
## TCGA.CV.6934.01 TCGA.CV.6934    Tumor
## TCGA.CV.6934.11 TCGA.CV.6934 Nontumor
## TCGA.CV.6935.01 TCGA.CV.6935    Tumor
## TCGA.CV.6935.11 TCGA.CV.6935 Nontumor
## TCGA.CV.6936.01 TCGA.CV.6936    Tumor
## TCGA.CV.6936.11 TCGA.CV.6936 Nontumor
## TCGA.CV.6938.01 TCGA.CV.6938    Tumor
## TCGA.CV.6938.11 TCGA.CV.6938 Nontumor

# 3.4.1 DEseq2的配对差异分析
dds_paired <- DESeqDataSetFromMatrix(countData = data,
                              colData = coldata_paired,
                              design = ~ subject + group)

## factor levels were dropped which had no samples

dds_paired$group <- relevel(dds_paired$group, ref = "Nontumor") # 指定哪一组作为对照

dds_paired <- DESeq(dds_paired)

## estimating size factors

## estimating dispersions

## gene-wise dispersion estimates

## mean-dispersion relationship

## final dispersion estimates

## fitting model and testing

nrDEG_paired <- as.data.frame(results(dds_paired))
nrDEG_paired = nrDEG_paired[order(nrDEG_paired$log2FoldChange),] # 按照logFC排序
nrDEG_paired[1:3,1:6]

##         baseMean log2FoldChange    lfcSE       stat       pvalue         padj
## PRH2   544288.00      -15.63870 2.509593  -6.231567 4.617911e-10 7.423411e-08
## STATH 1385580.80      -15.33676 1.292428 -11.866626 1.764414e-32 9.170836e-29
## MUC21   26826.61      -12.83806 1.239818 -10.354790 3.980606e-25 7.758699e-22

write.csv(nrDEG_paired,file = "./Output/TCGA_HNSCpaired_DESeq2_nrDEG.csv")

# 定义差异基因，用于后续的分析，例如火山图，GO和KEGG等
logFC_cutoff <- 2
nrDEG_paired$change <-  as.factor(ifelse(nrDEG_paired$padj < 0.05 & abs(nrDEG_paired$log2FoldChange) > logFC_cutoff,
                                  ifelse(nrDEG_paired$log2FoldChange > logFC_cutoff ,'UP','DOWN'),'NOT'))
table(nrDEG_paired$change)

## 
##  DOWN   NOT    UP 
##   753 14452   388

  # 绘制火山图
Volcano_data <- na.omit(nrDEG_paired)
Volcano_data$logP <- -log10(Volcano_data$padj)

Volcano_paired <- ggscatter(Volcano_data, x = "log2FoldChange", y = "logP",
                     color = "change",  
                     palette = c("#2f5688", "#BBBBBB", "#CC0000"),
                     size = 1,
                     font.label = 10, 
                     repel = T, 
                     xlab = "log2 FoldChange", 
                     ylab = "-log10 (pvalue)",   
                     title="Volcano Plot") + 
  theme_base() + 
  geom_hline(yintercept = 1.30, linetype="dashed") +
  geom_vline(xintercept = c(-logFC_cutoff,logFC_cutoff), linetype="dashed")+ 
  theme(legend.position = "right",plot.title = element_text(size = 14,color="black",hjust = 0.5))
Volcano_paired

ggsave("./Output/Volcano_paired.PDF",Volcano_paired,width=6,height=5)

上述就是配对DEseq2差异分析的一个结果。接下来，我们看一下如果使用普通DEseq2差异分析的结果，将会有哪些区别？

# 3.5 DEseq2的普通差异分析
dds_2 <- DESeqDataSetFromMatrix(countData = data,
                              colData = coldata_paired,
                              design = ~ group)
dds_2$condition <- relevel(dds_2$group, ref = "Nontumor") # 指定哪一组作为对照

dds_2 <- DESeq(dds_2)

## estimating size factors

## estimating dispersions

## gene-wise dispersion estimates

## mean-dispersion relationship

## final dispersion estimates

## fitting model and testing

nrDEG_2 <- as.data.frame(results(dds_2))
nrDEG_2 = nrDEG_2[order(nrDEG_2$log2FoldChange),] # 按照logFC排序
nrDEG_2[1:3,1:4]

##        baseMean log2FoldChange    lfcSE      stat
## PRH2   544288.0      -20.27705 2.722897 -7.446865
## STATH 1385580.8      -18.75144 2.069789 -9.059590
## PRH1   104842.6      -15.65276 1.820363 -8.598703

write.csv(nrDEG_2,file = "./Output/TCGA_HNSC_DESeq2_nrDEG.csv")

# 这里我还提取了标准化后的表达矩阵，可以用于后续的热图箱式图绘制等等，当然也可以使用FPKM或TPM代替
rld_paired <- vst(dds_paired)
Nor_expr_data_paired <- assay(rld_paired)
write.csv(Nor_expr_data_paired,file = "./Output/TCGA_HNSCpaired_Norexpr_data_paired.csv")

# 3.5 定义差异基因，用于后续的分析，例如GO和KEGG等
logFC_cutoff <- 2
nrDEG_2$change <-  as.factor(ifelse(nrDEG_2$padj < 0.05 & abs(nrDEG_2$log2FoldChange) > logFC_cutoff,
                                  ifelse(nrDEG_2$log2FoldChange > logFC_cutoff ,'UP','DOWN'),'NOT'))
table(nrDEG_2$change)

## 
##  DOWN   NOT    UP 
##   626 14093   375

# 3.6 绘制火山图
Volcano_data_2 = na.omit(nrDEG_2)
Volcano_data_2$logP <- -log10(Volcano_data_2$padj)
Volcano_2 <- ggscatter(Volcano_data_2, x = "log2FoldChange", y = "logP",
                     color = "change",  
                     palette = c("#2f5688", "#BBBBBB", "#CC0000"),
                     size = 1,
                     font.label = 10, 
                     repel = T, 
                     xlab = "log2 FoldChange", 
                     ylab = "-log10 (pvalue)",   
                     title="Volcano Plot") + 
  theme_base() + 
  geom_hline(yintercept = 1.30, linetype="dashed") +
  geom_vline(xintercept = c(-logFC_cutoff,logFC_cutoff), linetype="dashed")+ 
  theme(legend.position = "right",plot.title = element_text(size = 14,color="black",hjust = 0.5))
Volcano_2


ggsave("./Output/Volcano_2.PDF",Volcano_2 ,width=6,height=5)

可以看到配对的DEseq2差异分析有753基因下调，388基因上调；而常规的差异分析有626个基因下调，375基因上调，这个结果看似差不多 (生信宝典注：这里应该做个Venn图比较下的，除了数量，具体基因是否一致。可参考高通量数据中批次效应的鉴定和处理 - 系列总结和更新)。

实际上，笔者曾处理过一个项目，用配对和不用配对，得到的结果相差甚远，主要还是体现在火山图、KEGG和GO分析等需要输入差异基因的分析中。但好像做GSEA这样纳入全部基因的分析，并不会受到太大的影响，这点有机会再探讨。

生信宝典注 (具体见一文掌握GSEA，超详细教程！)：

  1. 如果用log2FC做GSEA结果是肯定会有差别的，影响多大不好说。

  2. 如果用标准化后的矩阵做GSEA结果没差别，因为配对分析没有影响到表达标准化，而只是在差异分析步骤起作用。

绘制单基因配对t检验箱式图

4.1 提取PRH2基因的表达量

## 4. PRH2表达量的配对箱式图绘制
PRH2_expression <- t(paired_expr_data["PRH2",] ) %>% data.frame(PRH2 = .)
PRH2_expression$group <- paired_phe$group
PRH2_expression$paired <- paired_phe$subject

PRH2_expression = PRH2_expression[order(as.numeric(paired_phe$subject)),] #这里需要按照配对信息排序好，以做后需的配对差异分析

4.2 开始绘图

  # 这个是我自己写的一个ggplot2的主题，可以自定义修改其中的参数
if(T){mytheme <- theme(plot.title = element_text(size = 12,color="black",hjust = 0.5),
                       axis.title = element_text(size = 12,color ="black"), 
                       axis.text = element_text(size= 12,color = "black"),
                       panel.grid.minor.y = element_blank(),
                       panel.grid.minor.x = element_blank(),
                       panel.grid=element_blank(),
                       legend.position = "none",
                       legend.text = element_text(size= 12),
                       legend.title= element_text(size= 12)
)} 

  # 配对t检验箱式图图
box_paired <- ggplot(PRH2_expression, aes(x = group,y = PRH2,fill= group)) +
  geom_boxplot(aes(fill = group),width= 1) +
  geom_point(size= 1, alpha=0.6) +
  geom_line(aes(group = paired), colour="black", linetype="11",size=0.5)+ 
  scale_fill_manual(values = c("#1CB4B8", "#EB7369"))+
  scale_y_continuous(limits=c(0, 25),breaks =seq(0,20,5),expand = c(0, 0) )+
  labs(y="log2 Expression",x=NULL,title = "TCGA")+ 
  theme_bw()+ mytheme +
  stat_compare_means(
    label = "p.format",
    method = "t.test",
    paired = T,
    size = 4, 
    label.x = 1.5,
    label.y = 22
)

box_paired

推荐一款高颜值免费在线SCI绘图工具~~~ 可以直接做配对箱线图

ggsave("./Output/TCGA_PRH2_paired_TumorvsNormal.PDF",
       box_paired ,height = 8 ,width = 6,units = "cm")

sessionInfo()

## R version 4.0.2 (2020-06-22)
## Platform: x86_64-w64-mingw32/x64 (64-bit)
## Running under: Windows 10 x64 (build 14393)
## 
## Matrix products: default
## 
## locale:
## [1] LC_COLLATE=Chinese (Simplified)_China.936  LC_CTYPE=Chinese (Simplified)_China.936   
## [3] LC_MONETARY=Chinese (Simplified)_China.936 LC_NUMERIC=C                              
## [5] LC_TIME=Chinese (Simplified)_China.936    
## 
## attached base packages:
## [1] parallel  stats4    stats     graphics  grDevices utils     datasets  methods   base     
## 
## other attached packages:
##  [1] ggthemes_4.2.0              ggpubr_0.4.0                ggplot2_3.3.2              
##  [4] dplyr_1.0.0                 DESeq2_1.28.1               SummarizedExperiment_1.18.2
##  [7] DelayedArray_0.14.0         matrixStats_0.56.0          Biobase_2.48.0             
## [10] GenomicRanges_1.40.0        GenomeInfoDb_1.24.2         IRanges_2.22.2             
## [13] S4Vectors_0.26.1            BiocGenerics_0.34.0        
## 
## loaded via a namespace (and not attached):
##  [1] tidyr_1.1.0            bit64_0.9-7            splines_4.0.2          carData_3.0-4         
##  [5] blob_1.2.1             cellranger_1.1.0       GenomeInfoDbData_1.2.3 yaml_2.2.1            
##  [9] pillar_1.4.6           RSQLite_2.2.0          backports_1.1.7        lattice_0.20-41       
## [13] glue_1.4.1             digest_0.6.25          RColorBrewer_1.1-2     ggsignif_0.6.0        
## [17] XVector_0.28.0         colorspace_1.4-1       htmltools_0.5.0        Matrix_1.2-18         
## [21] XML_3.99-0.4           pkgconfig_2.0.3        broom_0.7.0            haven_2.3.1           
## [25] genefilter_1.70.0      zlibbioc_1.34.0        purrr_0.3.4            xtable_1.8-4          
## [29] scales_1.1.1           openxlsx_4.1.5         rio_0.5.16             BiocParallel_1.22.0   
## [33] tibble_3.0.2           annotate_1.66.0        farver_2.0.3           generics_0.1.0        
## [37] car_3.0-8              ellipsis_0.3.1         withr_2.2.0            readxl_1.3.1          
## [41] survival_3.1-12        magrittr_1.5           crayon_1.3.4           memoise_1.1.0         
## [45] evaluate_0.14          rstatix_0.6.0          forcats_0.5.0          foreign_0.8-80        
## [49] tools_4.0.2            data.table_1.12.8      hms_0.5.3              lifecycle_0.2.0       
## [53] stringr_1.4.0          munsell_0.5.0          locfit_1.5-9.4         zip_2.1.1             
## [57] AnnotationDbi_1.50.1   compiler_4.0.2         rlang_0.4.6            grid_4.0.2            
## [61] RCurl_1.98-1.2         labeling_0.3           bitops_1.0-6           rmarkdown_2.3         
## [65] gtable_0.3.0           abind_1.4-5            DBI_1.1.0              curl_4.3              
## [69] R6_2.4.1               knitr_1.29             bit_1.1-15.2           stringi_1.4.6         
## [73] Rcpp_1.0.5             vctrs_0.3.1            geneplotter_1.66.0     tidyselect_1.1.0      
## [77] xfun_0.15

参考资料

  1. 《利用R代码从UCSC XENA下载mRNA, lncRNA, miRNA表达数据并匹配临床信息》

  2. 《差异分析|DESeq2完成配对样本的差异分析》

  3. 《R数据科学》

  4. 数据和代码下载：

    https://gitee.com/ct5869/shengxin-baodian/tree/master/TCGA

  5. 作者：赵法明

    编辑：生信宝典

往期精品(点击图片直达文字对应教程)

机器学习

后台回复“生信宝典福利第一波”或点击阅读原文获取教程合集

生信宝典
关注
  • 14
    点赞
  • 踩
  • 156
    收藏
    觉得还不错? 一键收藏
  • 打赏
    打赏
  • 6
    评论
R语言_差异表达分析
10-21
介绍R语言用于统计分析差异表达，来源课程内部讲义，注重实战！
生信分析论文套路R语言代码
12-05
TCGA GEO数据处理，基因注释，差异分析，GO，KEGG，GSEA富集分析，肿瘤突变负荷，免疫浸润，LASSO回归，随机森林，SVM-RFE，COX回归，WGCNA网络，共识聚类分析，药物敏感性分析，干性指数，免疫浸润指数，预后模型等维度的相关R语言代码
TCGA数据分析-差异表达.R.r
11-30
TCGA数据分析-差异表达.R.r
TCGA数据分析.mp4
06-18
TCGA数据分析
生物信息学概论——聚类分析TCGA-BRCA数据.zip
07-05
利用R软件或其他数据分析语言，进行该数据的聚类分析。利用层次聚类，对该组数据样本按照基因表达水平进行聚类，看聚类效果如何。即是否能够按照基因表达水平，将病人进行分类。距离可以选择average。注，R中有相应的聚类函数，请利用并尽可能输出图示（如heatmap），表明你的结果。实现PCA，并利用你实现的PCA对该组数据的基因表达进行降维处理。请选择你认为合适的主成分数目，给出原因，再次对病人依据你给的特征进行聚类，并与1比较。 数据文件说明:GeneMatrix.txt：基因表达值文件，含有行名和列名，一行为一个基因，一列为一个病人clinical_data.txt：记录了病人的若干信息，每一行为一个病人，病人的编号和GeneMatrix.txt中的相同。GeneMatrix中病人只涵盖了这里的一部分，注意，在病人的若干描述中，有一项为ER_Status_nature2012，可以根据这个对病人进行分类，你可以按照这个分类标准，对你的聚类进行一定的评估，看结果是否符合预期。 参考：https://blog.csdn.net/newlw/article/details/124661950
R语言零基础基因/数据差异分析（一）
Frms的博客
07-15 1万+
文章目录介绍环境搭建软件下载结果展示基因数据下载流程基因数据处理利用GEO分析绘制拟火山图 注意，本 系列 有连贯性，每一步都很详细，每一步都很重要，请耐心读完！！ 介绍 本系列文主要依据真实论文制图流程，详细说明制图过程， 其中包括： 1. 基因数据下载 2. 制图所需数据格式 3. 火山图制作流程 4. 聚类热图制作流程 环境搭建 软件下载 移步至此学习 结果展示 基因数据处理 注意删除末行注释 基因数据下载流程 以GSE137578基因为例，先下载 如图所示文件，并解压如图。
R语言limma包差异表达分析
m0_64111502的博客
11-11 7320
本文主要是对差异表达分析进行演示，其中利用的是R包limma.
R语言——基于DESeq2做基因表达差异分析
codyliu_bo的博客
10-28 4228
本记录主要是使用DESeq2 包进行差异包表达分析。R 语言是4.2.0版本，对应的BiocManager是3.15版本。我们采用的是BiocManager安装DESeq2包。检查数据Counts文件与coldata数据是否匹配。注意部分依赖包是由R4.2.1建造，可能会有问题。过滤在所有重复样本中小于1的基因。
R语言零基础基因/数据差异分析（二）
Frms的博客
07-15 6203
文章目录结果展示安装ggplot2包制作方法 注意，本 系列 有连贯性，每一步都很详细，每一步都很重要，请耐心读完！！ 结果展示 安装ggplot2包 如图操作 找到并勾上即可。 制作方法 关于 基因啊选定标准，即logFC和FDR，请仔细阅读代码修改即可，这里给出的标准是： |log2(FC)| > 1 且 FDR < 0.01 ```r #选择文件 df=read.csv(file.choose(), header = T #是否有标题，T表示有，F反之 ) #加载包（反正多余
差异基因分析R代码
01-06
资源有用，希望能够对大家有帮助。好资源大家gongxiang
TCGA数据下载及全流程分析（更新中）
01-20
一、GDCRNATools包下载 首先下载R包 if (!requireNamespace(BiocManager, quietly = TRUE)) install.packages(BiocManager) BiocManager::install(GDCRNATools) 下载好之后运行，以COAD数据为例，下载： library(GDCRNATools) gdcRNADownload(project.id = TCGA-COADREAD, data.type = RNAseq, write.ma
新版TCGA数据整理，通过R语言实现TCGA数据整理
04-05
新版TCGA数据的整理，将TCGA单个转录组数据整理成列为患者，行为基因的矩阵。可以通过修改参数整理Count、FPKM、TPM数据，根据自己的需要进行选择。
R语言绘制圈图、环形热图可视化基因组实战：展示基因数据比较
拓端研究室TRL
09-27 293
原文链接：http://tecdat.cn/?p=23891可以使用环状图形展示基因数据比较。可以添加多种图展信息，如热图、散点图等。本文目标:可视化基因组数据制作环形热图环形热图很漂亮。可以通过R来实现环形热图。首先，让我们生成一个随机矩阵，并将其随机分成五组。mat1=rbind(cbind(matrix(rnorm(50*5,mean=1),nr=50), ...
R语言差异检验：t检验
热门推荐
tidyss
11-09 1万+
文章目录@[toc]单样本t检验适用条件具体计算公式R语言示例独立样本t检验适用条件具体计算公式R语言示例配对样本t检验适用条件具体计算公式R语言示例 t检验(student t检验)是应用t分布的特征，将t作为检验的统计量来进行统计推断方法。它对样本要求较小(例如n＜30)。 主要用途： 样本均数与总体均数的差异比较 两样本均数的差异比较 分类： 单样本t检验 独立样本t检验 配对样本t...
TCGA数据库学习二：差异基因分析
weixin_73362123的博客
03-31 3991
在分析TCGA阿虎句酷里的RNA-seq数据之前，先了解TCGA样本的id名称。TCGA样本id中第14-15位代表分组信息，01-09是tumor(肿瘤)，10-29是normal。
百日筑基篇——差异基因分析DESeq2包（R语言初识六）
2301_78630677的博客
07-14 1540
今天，整理一下如何使用R语言进行基因差异分析。主要会讲述有关的几个扩展包以及可视化，希望以这种方式巩固所学。基因差异分析的学习，今天就到这里，理解这背后的统计学原理是十分必要的，同时要一开始就确定想要通过基因差异分析，分析出什么，目标明确。 醉里挑灯看剑，梦回吹角连营。 -2023-7-8 筑基篇。
百日筑基篇——差异基因分析Limma包（R语言初识七）
2301_78630677的博客
07-15 4217
Limma包是针对微阵列数据和RNA-seq数据开发的一个强大的R软件包，其提供了一套统计模型和工具来鉴定差异表达基因，是进行差异基因分析的常用方法之一。 世上无难事，只要肯攀登 -2023-7-9 筑基篇
【 书生·浦语大模型实战营】作业（六）：Lagent & AgentLego 智能体应用搭建
最新发布
专注大数据与人工智能技术分享，欢迎私信加群互相学习！
05-03 1370
本篇笔记内容主要为 【书生·浦语大模型实战营】作业（六）：Lagent & AgentLego 智能体应用搭建智能体的作业 ，主要对智能体应用理论部分及具体搭建流程详细介绍，并分别进行实战应用及可视化展示，欢迎大家交流学习！
SpringBoot中实现发送邮件
hac1322的博客
04-30 807
当你添加了spring-boot-starter-mail依赖后，Spring Boot会自动配置JavaMailSender实例，并根据application.yml文件中的属性来配置这个实例。你可以直接在需要发送邮件的地方通过@Autowired注解将JavaMailSender实例注入到你的类中，然后使用它来发送邮件。：首先，需要在你的pom.xml文件中添加Spring Boot的邮件发送器依赖。，简化了在Spring Boot应用程序中发送电子邮件的设置过程。Spring Boot的。
R语言 差异基因分析
08-30
在R语言中进行差异基因分析一般可以使用DEseq2包。首先，需要对原始数据进行预处理，包括去除低表达的基因和样本间的批次效应。然后，可以使用DEseq2进行配对差异分析。在配对差异分析中，保留至少在75%的样本中都有表达的基因。 DEseq2的配对差异分析结果可以显示下调和上调的基因数量。例如，可以发现有753个基因下调和388个基因上调。与常规的差异分析相比，配对差异分析的基因数量可能类似，但具体的基因是否一致需要通过Venn图进行比较。 值得注意的是，如果使用标准化后的矩阵进行GSEA分析，配对分析不会影响表达标准化，而只会在差异分析步骤起作用。此外，还可以绘制单基因配对t检验的箱式图来展示差异基因的表达情况。<span class="em">1</span><span class="em">2</span><span class="em">3</span> #### 引用[.reference_title] - *1* *2* *3* [R语言进行TCGA配对样本差异基因分析](https://blog.csdn.net/qazplm12_3/article/details/115339797)[target="_blank" data-report-click={"spm":"1018.2226.3001.9630","extra":{"utm_source":"vip_chatgpt_common_search_pc_result","utm_medium":"distribute.pc_search_result.none-task-cask-2~all~insert_cask~default-1-null.142^v93^chatsearchT3_2"}}] [.reference_item style="max-width: 100%"] [ .reference_list ]

“相关推荐”对你有帮助么？

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论 6
添加红包

请填写红包祝福语或标题

个

红包个数最小为10个

元

红包金额最低5元

当前余额3.43元 前往充值 >
需支付：10.00元
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包

打赏作者

生信宝典

你的鼓励将是我创作的最大动力

¥1 ¥2 ¥4 ¥6 ¥10 ¥20
扫码支付：¥1
获取中
扫码支付

您的余额不足，请更换扫码支付或充值

打赏作者

实付元
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明：

1.余额是钱包充值的虚拟货币，按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载，可以购买VIP、付费专栏及课程。

余额充值