专栏首页生信小驿站多分组差异分析解决方案(1)循环T检验

多分组差异分析解决方案(1)循环T检验

主要方法:将其中某一组设置为实验组,其余几组统一设置为对照组。

  • 第一步读取数据,合并表达矩阵和分组文件
#===========================================================================


#===========================================================================


rm(list = ls(all.names = TRUE))

options(stringsAsFactors = F)

library(Matrix)



setwd('D:\\SCIwork\\F38KRT\\s2')


data <- read.csv('cdata.csv', header = T, row.names = 1)

data <- as.data.frame(t(data))

data[1:4,1:4]


normalized<-function(y) {
  
  x<-y[!is.na(y)]
  
  x<-(x - min(x)) / (max(x) - min(x))
  
  y[!is.na(y)]<-x
  
  return(y)
}

db <-  as.data.frame(apply(data,2,normalized))


data <- db

data$sample <- rownames(data)

data$sample <- chartr(old='.', new='-', x=data$sample)


setwd('D:\\SCIwork\\F38KRT\\s3')

group <- read.csv('group2.csv', header = T)

names(group)[1] <- 'sample'

group$sample <- chartr(old='.', new='-', x=group$sample)

group <- subset(group, select=c("sample", "group"))

group$subtype <- group$group

group$group  <- NULL




dt <- merge(group, data, by='sample')

dt[1:4,1:4]

dt$sample <- NULL

table(dt$subtype)

dt_total <- dt

输入文件cdata经过T转置为data后如下所示:

上面的代码包括了,对表达量归一化的代码。

normalized<-function(y) { x<-y[!is.na(y)] x<-(x - min(x)) / (max(x) - min(x)) y[!is.na(y)]<-x return(y)}

经过这个代码,样本的表达量将会被转化到0-1之间的数值。

  • 将subtype1设置为exp组,其余两组(subtype2,和subtype3)设置为con组。
#===========================================================================


#===========================================================================


dt$subtype <- ifelse(dt$subtype == 'Subtype1', 'Exp', 'Con')

table(dt$subtype)

dt <- dt[order(dt$subtype), ]

dt[1:4,1:4]

dt_Con <- subset(dt, dt$subtype == 'Con')

dt_Con[1:4,1:4]

dt_Exp <- subset(dt, dt$subtype == 'Exp')

dt_Exp[1:4,1:4]




dt_Con$subtype <- paste0(dt_Con$subtype, rownames(dt_Con))

rownames(dt_Con) <- dt_Con$subtype

dt_Con$subtype <- NULL

dt_Con <- as.data.frame(t(dt_Con))




dt_Exp$subtype <- paste0(dt_Exp$subtype, rownames(dt_Exp))

rownames(dt_Exp) <- dt_Exp$subtype

dt_Exp$subtype <- NULL

dt_Exp <- as.data.frame(t(dt_Exp))


Pvalue<-c(rep(0,nrow(dt_Con)))

log2_FC<-c(rep(0,nrow(dt_Con)))

for(i in 1:nrow(dt_Con)){
  
  y=t.test(as.numeric(dt_Con[i,]),as.numeric(dt_Exp[i,]))
  Pvalue[i] <- y$p.value
  log2_FC[i] <-log2(mean(as.numeric(dt_Exp[i,]))/(mean(as.numeric(dt_Con[i,]))))
  
}


library(dplyr)

library(tidyr)

library(tibble)

# 对p value进行FDR校正
fdr=p.adjust(Pvalue, "BH") 
# 在原文件后面加入log2FC,p value和FDR,共3列;
out<- as.data.frame(cbind(log2_FC,Pvalue,fdr))
out$gene <- rownames(dt_Con)
# out <- out %>%
#   dplyr::filter(log2_FC > 0.5 & Pvalue < 0.05)


setwd('D:\\SCIwork\\F38KRT\\s5')

write.csv(out, file = 'out_S1.csv')

循环T检验后求取差异基因的差异倍数和P值。

  • 同样的逻辑,分别求取subtype2和subtype3的差异基因
#===========================================================================


#===========================================================================



dt <- dt_total

dt$subtype <- ifelse(dt$subtype == 'Subtype2', 'Exp', 'Con')

table(dt$subtype)

dt <- dt[order(dt$subtype), ]

dt[1:4,1:4]

dt_Con <- subset(dt, dt$subtype == 'Con')

dt_Con[1:4,1:4]

dt_Exp <- subset(dt, dt$subtype == 'Exp')

dt_Exp[1:4,1:4]




dt_Con$subtype <- paste0(dt_Con$subtype, rownames(dt_Con))

rownames(dt_Con) <- dt_Con$subtype

dt_Con$subtype <- NULL

dt_Con <- as.data.frame(t(dt_Con))




dt_Exp$subtype <- paste0(dt_Exp$subtype, rownames(dt_Exp))

rownames(dt_Exp) <- dt_Exp$subtype

dt_Exp$subtype <- NULL

dt_Exp <- as.data.frame(t(dt_Exp))


Pvalue<-c(rep(0,nrow(dt_Con)))

log2_FC<-c(rep(0,nrow(dt_Con)))

for(i in 1:nrow(dt_Con)){
  
  y=t.test(as.numeric(dt_Con[i,]),as.numeric(dt_Exp[i,]))
  Pvalue[i] <- y$p.value
  log2_FC[i] <-log2(mean(as.numeric(dt_Exp[i,]))/(mean(as.numeric(dt_Con[i,]))))
  
}



# 对p value进行FDR校正
fdr=p.adjust(Pvalue, "BH") 
# 在原文件后面加入log2FC,p value和FDR,共3列;
out<- as.data.frame(cbind(log2_FC,Pvalue,fdr))
out$gene <- rownames(dt_Con)



setwd('D:\\SCIwork\\F38KRT\\s5')

write.csv(out, file = 'out_S2.csv')







#===========================================================================


#===========================================================================



dt <- dt_total

dt$subtype <- ifelse(dt$subtype == 'Subtype3', 'Exp', 'Con')

table(dt$subtype)

dt <- dt[order(dt$subtype), ]

dt[1:4,1:4]

dt_Con <- subset(dt, dt$subtype == 'Con')

dt_Con[1:4,1:4]

dt_Exp <- subset(dt, dt$subtype == 'Exp')

dt_Exp[1:4,1:4]




dt_Con$subtype <- paste0(dt_Con$subtype, rownames(dt_Con))

rownames(dt_Con) <- dt_Con$subtype

dt_Con$subtype <- NULL

dt_Con <- as.data.frame(t(dt_Con))




dt_Exp$subtype <- paste0(dt_Exp$subtype, rownames(dt_Exp))

rownames(dt_Exp) <- dt_Exp$subtype

dt_Exp$subtype <- NULL

dt_Exp <- as.data.frame(t(dt_Exp))


Pvalue<-c(rep(0,nrow(dt_Con)))

log2_FC<-c(rep(0,nrow(dt_Con)))

for(i in 1:nrow(dt_Con)){
  
  y=t.test(as.numeric(dt_Con[i,]),as.numeric(dt_Exp[i,]))
  Pvalue[i] <- y$p.value
  log2_FC[i] <-log2(mean(as.numeric(dt_Exp[i,]))/(mean(as.numeric(dt_Con[i,]))))
  
}



# 对p value进行FDR校正
fdr=p.adjust(Pvalue, "BH") 
# 在原文件后面加入log2FC,p value和FDR,共3列;
out<- as.data.frame(cbind(log2_FC,Pvalue,fdr))
out$gene <- rownames(dt_Con)



setwd('D:\\SCIwork\\F38KRT\\s5')

write.csv(out, file = 'out_S3.csv')

本文参与腾讯云自媒体分享计划,欢迎正在阅读的你也加入,一起分享。

我来说两句

0 条评论
登录 后参与评论

相关文章

  • 多分组差异分析解决方案(2)分批次差异基因后取交集

    主要方法:如果不同分组代表着一定的趋势,例如group1,group2,group3的样本严重程度越来越重。那么就可以求group1和group2的差异基因,g...

    用户1359560
  • 「R」R检验中的“数据是恆量”问题

    这是一般做基因差异表达分析在使用t检验或者其他统计检验中常出现的一个问题。之前我学习和自己分析时就遇到过,尝试使用判断的方式事先检查它是不是数据存在问题(这类数...

    王诗翔呀
  • elife: 写作及审稿中常见的十个统计错误

    本文列出了文献中出现的一些最常见的统计错误。这些错误的根源在于无效的实验设计、不恰当的分析或有缺陷的推理。作者对如何识别和解决这些错误为研究者和审稿人提供了建议...

    Listenlii-生物信息知识分享
  • 使用student’s T检验的未必是学生

    一直想整理一下统计方法在网站分析中的应用,刚好前几天遇到类似的问题,借这个机会整理一下网站分析中T检验的思路。在统计面前我们并没有生产方法,我们只是方法的搬...

    数据森麟
  • 「Workshop」第一期:我理解的(生信)数据分析核心基础

    我在简书和公众号上已经分享了很多之前学习的数据分析笔记和文章,覆盖了各方面的内容,数据分析方面以后不会再个人分享特别基础的东西了。接下来我会让师弟师妹们定期分享...

    王诗翔呀
  • 美团配送A/B评估体系建设与实践

    2019年5月6日,美团正式推出新品牌“美团配送”,发布了美团配送新愿景:“每天完成一亿次值得信赖的配送服务,成为不可或缺的生活基础设施。”现在,美团配送已经服...

    美团技术团队
  • 线性回归(二)-违背基本假设的情况和处理方法

    由线性回归(一)^1,我们通过数学中的极值原理推导出了一元线性回归的参数估计和多元线性回归的参数估计的拟合方程计算方法。同时为了检验拟合质量,我们引入了两种主要...

    EatRice
  • 干货 | 一个数据分析师眼中的数据预测与监控

    束开亮,携程大市场部BI团队,负责数据分析与挖掘。同济应用数学硕士,金融数学方向,法国统计学工程师,主修风险管理与金融工程。

    携程技术
  • 方差分析简介(结合COVID-19案例)

    我们正在应对一场空前规模的流行病。全世界的研究人员都在疯狂地试图开发一种疫苗或COVID-19的治疗方法,而医生们正试图阻止这种流行病席卷整个世界。

    磐创AI

扫码关注云+社区

领取腾讯云代金券