1. 相关性散点图展示
d3fdd706c6cf8c2f4ea2d390ecb3bbc.png
2. 代码块
options(warn = -1)
RNA = read.table('./细胞A/RNA.txt',header =T,sep = '\t')
RNA =RNA[,-4]
RNA1 <- RNA[!duplicated(RNA$gene_symbol),]
rownames(RNA1) <- RNA1[,1]
RNA1 <- RNA1[,-1]
RNA2 = RNA1[,c(4:7,1:3)]
RNA3 <- log(RNA2+1)
#dat1 = RNA2
#dat1<-dat1[!apply(dat1,1,function(x){sum(floor(x)==0)>1}),]
library(limma)
group <- c(rep("young",4),rep("aged",3))
group <- factor(group)
design <- model.matrix(~0 + group)
colnames(design) <- levels(group)
contrast.matrix <- makeContrasts(aged - young,
levels=design)
contrast.matrix
fit <- lmFit(RNA3,design)
fit2 <- contrasts.fit(fit, contrast.matrix)
fit2 <- eBayes(fit2)
allDiff1=topTable(fit2,adjust='fdr',coef=1,number=Inf)
allDiff = na.omit(allDiff1)
write.csv(allDiff, file = "young_vs_aged_rna.csv" )
allDiff = allDiff %>% tibble::rownames_to_column('id')
table(allDiff$adj.P.Val> 0.05)
### 分泌蛋白差异分析
pro = read.table('./细胞A/sepro.txt',header =T,sep = '\t')
pro1 <- pro[!duplicated(pro$Gene.name),]
rownames(pro1) <- pro1[,1]
pro1 <- pro1[,-1]
pro1 = pro1[,c(4:6,1:3)]
head(pro1)
pro1 <- log(pro1+1)
head(pro1)
pro1 = na.omit(pro1)
library(limma)
group <- c(rep("young",3),rep("aged",3))
group <- factor(group)
design <- model.matrix(~0 + group)
colnames(design) <- levels(group)
contrast.matrix <- makeContrasts(aged - young,
levels=design)
contrast.matrix
fit <- lmFit(pro1,design)
fit2 <- contrasts.fit(fit, contrast.matrix)
fit2 <- eBayes(fit2)
allDiff1=topTable(fit2,adjust='fdr',coef=1,number=Inf)
limma.na <- na.omit(allDiff1)
write.csv(limma.na, file = "young_vs_aged_pro.csv" )
library(tibble)
limma.na = limma.na %>% tibble::rownames_to_column('id')
head(limma.na,2)
combine= merge(allDiff,limma.na,
by.x="id",
by.y="id",
suffixes = c("_RNA","_Protein") ,
all.x=FALSE,
all.y=FALSE)
dim(combine)
write.csv(combine,"RNA_protein3.csv",row.names=FALSE)
library(dplyr)
library(ggplot2)
library(ggrepel)
colnames(combine)
data <- data.frame(combine[c(1,2,6,8,12)])
tail(data)
data1 = data %>%
mutate(part = case_when(
logFC_RNA > 0 & adj.P.Val_RNA < 0.05 & adj.P.Val_Protein < 0.05 & logFC_Protein > 0 ~ "mRNA and Protein enriched in Aged",
logFC_RNA < 0 & logFC_Protein < 0 & adj.P.Val_RNA < 0.05 & adj.P.Val_Protein < 0.05 ~ "mRNA and Protein enriched in Young",
adj.P.Val_Protein > 0.05 | data$adj.P.Val_RNA > 0.05 ~ ' No sig(p.adj > 0.05)' ,
TRUE ~ 'mRNA and protein with opposite exprssion '))
table(data1$part)
mycolor <- c("#8B0000","#191970","black","gray80")
cor = round(cor(data$logFC_RNA,data$logFC_Protein,method = 'spearman'),2)
lab = paste('Spearman','\n',"correlation=",cor,sep="", '\n', 'p-value < 2.2e-16')
cor.test(data$logFC_RNA,data$logFC_Protein,method = 'spearman')
data1$part = factor(data1$part,levels= c("mRNA and Protein enriched in Aged",
"mRNA and Protein enriched in Young",
'mRNA and protein with opposite exprssion ',
' No sig(p.adj > 0.05)'))
p1 = ggplot(data1,aes(logFC_Protein,logFC_RNA,color=part))+
geom_point(size=2)+
scale_colour_manual(name="",values=alpha(mycolor,0.7))+
geom_hline(yintercept = c(0),
size = 0.5,
color = "grey40",
lty = "solid")+
geom_vline(xintercept = c(0),
size = 0.5,
color = "grey40",
lty = "solid")
reg<-lm(formula = logFC_Protein ~ logFC_RNA,
data=data1)
#get intercept and slope value
coeff<-coefficients(reg)
intercept<-coeff[1]
slope<- coeff[2]
p1
library(dplyr) # 用于数据处理
top_5 <- bind_rows(
data1 %>%
filter(part == 'mRNA and Protein enriched in Aged') %>%
arrange(adj.P.Val_RNA, desc(logFC_RNA)) %>%
head(5),
data1 %>%
filter(part == 'mRNA and Protein enriched in Young') %>%
arrange(adj.P.Val_RNA, desc(logFC_RNA)) %>%
head(5)
)
diy = data1 %>% filter(id %in% c('Vcan'))
plog = rbind(top_5,diy)
p1+
geom_label_repel(data = plog,
aes(label = id),
box.padding = unit(0.5, "lines"),
point.padding = unit(0.8, "lines"), segment.color = "black", show.legend = FALSE) +
theme_bw()+
theme(legend.position = c(.3,.9)) +
xlab('Cell Proteome Level\nLog2Fold Change(Aged/Young)')+
ylab('Cell RNA Level\nLog2Fold Change(Aged/Young)')+
theme(text = element_text(size = 13)) +
#ggtitle("Log2 fold change: Young vs Aged") +
theme(plot.title = element_text(hjust = 0.5))+
geom_abline(intercept = intercept, slope = slope, color="grey",
linetype="solid", size=0.5)+
theme(legend.background=element_rect(fill=rgb(1,1,1,alpha=0.001),colour=NA))+
annotate("text", x=2, y=-3, label=lab, color = "black")
#geom_text(x = 2, y = -3, label = lab,color = 'black',family = 'sans')
ggsave('A-cell.pdf')
