Month: March 2021

data <- read.delim(‘mouse_gc_log2.xls’, row.names = 1, sep = ‘\t’, check.names = FALSE)#行为基因ID列为样本的表达矩阵 data[order(rowSums(data),decreasing=T)[1:3000],] #筛选表达量高的前3000个基因 data[rowSums(data)>1,] #过滤掉表达量低的基因 data1=rowMeans(data)>1 #按平均数筛选 data2=rowSums(data>0)>6 #表达量不为0的样品个数筛选 data=data[data1 & data2,] #联合一下 data$cv <- apply(data, 1, function(x){ sd(x)/mean(x)*100 }) data_df <- data[order(data$cv, decreasing = T)[1:3000], 1:15]#筛选变异系数最大的3000个基因

library(edgeR) library(ggplot2) data <- read.delim(‘../mouse_TPM_notCross.matrix.xls’, row.names = 1, sep = ‘\t’, check.names = FALSE) group <- factor(rep(c(‘A’, ‘B’, ‘C’), each = 3)) y <- DGEList(counts = data, group = group) y <- calcNormFactors(y) logcpm <- cpm(y, prior.count=3, log=TRUE) write.table(logcpm, file=”mouse_TMM.xls”, sep=”\t”, quote=F, row.names=T, col.names=T) #dgList <- estimateCommonDisp(y) #dgList <- estimateTagwiseDisp(dgList) #norm_counts.table <- t(t(dgList$pseudo.counts)*(dgList$samples$norm.factors)) #write.table(norm_counts.table, file=”mouse_gc_pareto_TMM.xls”,…

什么是批次效应（batch effect）？不同平台的数据，同一平台的不同时期的数据，同一个样品不同试剂的数据，以及同一个样品不同时间的数据等等都会产生一种batch effect 。这种影响如果广泛存在应该被足够重视，否则会导致整个实验和最终的结论失败。比对实验组和对照组，不同的处理是患病和不患病（测序时，先测得疾病，然后测得正常），然后你通过分析，得到很多差异表达的基因。现在问题来了，这个差异表达的结果是和你要研究的因素有关，还是时间有关，这个问题里时间就会成为干扰实验结果的因素，这个效应就是batch effect。 library(limma) data <- read.table(“mouse_gc.xls”,sep=”\t”,header = TRUE) batch1 <- c(rep(‘fish’,9),rep(‘mouse’,9)) batch1 <- as.factor(batch1) design <- model.matrix(~0 + batch1)#data=normalizeBetweenArrays(data)#如果组内的中位数不在同一条水平线上，现在该参数校正，再使用下面的批次校正 new_data <- removeBatchEffect(data[,2:19], batch = batch1) boxplot(data[,2:19]) boxplot(new_data) write.table(new_data, file=”removeBatch.xls”, sep=”\t”, quote=F, row.names=T, col.names=T) limma结果出现负值请参考：https://cloud.tencent.com/developer/article/1680305参考链接：https://www.omicsclass.com/article/1113

library(“preprocessCore”) library(“ggplot2”) data <- read.table(“mouse_gc.xls”,sep=”\t”,header = TRUE) rownames(data) <- data$GeneID #rownames(data) <- data[,1] #data_mat <- data.matrix(data[,-1]) quantile_normalisation <- function(df){ df_rank <- apply(df,2,rank,ties.method=”min”) df_sorted <- data.frame(apply(df, 2, sort)) df_mean <- apply(df_sorted, 1, mean) index_to_mean <- function(my_index, my_mean){ return(my_mean[my_index]) } df_final <- apply(df_rank, 2, index_to_mean, my_mean=df_mean) rownames(df_final) <- rownames(df) return(df_final) } new_data <- quantile_normalisation(data[,2:19]) boxplot(data) boxplot(new_data) write.table(new_data,…