-
Notifications
You must be signed in to change notification settings - Fork 2
/
Copy path01_ABCD.R
69 lines (57 loc) · 2.6 KB
/
01_ABCD.R
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
rm(list=ls());
options(stringsAsFactors = F);
## load data
load('../data/cd4_dat.Rdata')
load('../data/cd8_dat.Rdata')
load('../data/NK_dat.Rdata')
load('../data/other_dat.Rdata')
load('../data/b_dat.Rdata')
source('supp_func.R');
## merge data
dat <- dplyr::bind_cols(cd4_dat,cd8_dat,b_dat,NK_dat,other_dat)
## gene annotation
ann <- read.table("../data/01/gene4matrix_Seurat.txt", header=T, row.names=1);
sum(row.names(dat)==row.names(ann))==nrow(ann);
## relationships of 3 pathway genes of KEGG
gene <- read.table("../data/01/hsa04662_f.txt", header=T);
g.u <- unique(c(gene[,1], gene[,2]));
dat.i <- which(ann[,1] %in% g.u);
dat.s <- dat[dat.i,];
## count zero number of genes among cells
z1 <- apply(dat, 1, function(x) sum(x==0));
z1 <- z1/ncol(dat) * 100;
z2 <- apply(dat.s, 1, function(x) sum(x==0));
z2 <- z2/ncol(dat.s) * 100;
z3 <- apply(dat, 2, function(x) sum(x==0));
z3 <- z3/nrow(dat) * 100;
z4 <- apply(dat.s, 2, function(x) sum(x==0));
z4 <- z4/nrow(dat.s) * 100;
pdf("1_ABCD.pdf", 23, 20);
layout(matrix(c(1,2,3,4), nrow = 2, ncol = 2, byrow = TRUE));
par("mar"=c(7, 7, 7, 7));
cols <- rainbow(10, alpha=0.8);
hist(z1, col=cols[7], border=cols[7], breaks=100, xlim=c(70, 100),
font.main = 1,main="21,430 full genes; 15,973 single cells", xlab="",
ylab = '',cex.main=3,cex.lab=3.5, cex.axis = 2.5);
title(ylab="Number of Genes", line=4, cex.lab=3.5, family = "sans")
title(xlab="% of Zero Value", line=4.5, cex.lab=3.5, family = "sans")
fig_label('A', pos='topleft',cex=5)
hist(z3, col=cols[7], border=cols[7], breaks=100, xlim=c(70, 100),
font.main = 1,main="15,973 single cells; 21,430 genes", xlab="",
ylab = '',cex.main=3,cex.lab=3.5, cex.axis= 2.5);
title(ylab="Number of Cells", line=4, cex.lab=3.5, family = "sans")
title(xlab="% of Zero Value", line=4.5, cex.lab=3.5, family = "sans")
fig_label('B', pos='topleft',cex=5)
hist(z2, col=cols[7], border=cols[7], breaks=100, xlim=c(70, 100),
font.main = 1,main="347 pathway genes with 15,973 single cells", xlab="",ylab = '',
cex.main=3,cex.lab=3.5, cex.axis= 2.5);
title(ylab="Number of Genes", line=4, cex.lab=3.5, family = "sans")
title(xlab="% of Zero Value", line=4.5, cex.lab=3.5, family = "sans")
fig_label('C', pos='topleft',cex=5)
hist(z4, col=cols[7], border=cols[7], breaks=100, xlim=c(70, 100),
font.main = 1,main="15,973 single cells with 347 pathway genes", xlab="",ylab = '',
cex.main=3,cex.lab=3.5, cex.axis= 2.5);
title(ylab="Number of Cells", line=4, cex.lab=3.5, family = "sans")
title(xlab="% of Zero Value", line=4.5, cex.lab=3.5, family = "sans")
fig_label('D', pos='topleft',cex=5)
dev.off();