-
Notifications
You must be signed in to change notification settings - Fork 1
/
simplecheck.Rmd
174 lines (154 loc) · 8.13 KB
/
simplecheck.Rmd
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
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
---
title: "Laterally confined growth for stem cell engineering: Looking at the trascriptome"
author: "Saradha Venkatachapathy"
output:
html_document:
fig_caption: yes
fig_height: 6
highlight: tango
toc: yes
---
# Introduction
In the paper, we present a novel assay for reprogramming fibroblasts into iPSC-like cells in the absence of exogeous factors. Briefly, mouse embryonic fibroblasts (NIH3T3) were cultured such that their growth was confined on rectangular islands. After 10 days of such latterally confined growth, we observed that they acquired stem cell characteristics. We perfored RNASeq and obtained the transciption profile of cells after 3 hours(S1), 3 days(S2), 6 days(S3) and 10 days(S4) of laterally confined growth on micropatterns. Here we present proof that latterally confined growth induces nuclear reprogramming.
#Characterizing the expression of candidate genes
We read in the gene expression counts and check the expression of candidate genes. Note, the candidate gene and the important gene lists were complied from Literature Review.
Labels: Mesenchymal (Me), Embryonic Stem cell(ES) and induced Pluripotent Stem Cells(iPSC).
```{r Read in the gene counts, message=FALSE, warning=FALSE}
rnaseq_mapped_genes<-read.csv("./allcombined_gene_level_expression_mapped_to_gene_name.csv")
rnaseq_mapped_genes<-rnaseq_mapped_genes[,-1]
rnaseq_mapped_genes[,79]<-apply(rnaseq_mapped_genes[,75:78],1,mean)
rnaseq_mapped_genes[,80]<-apply(rnaseq_mapped_genes[,75:78],1,sd)
colnames(rnaseq_mapped_genes)[79:80]<-c("Pop_Mean","Pop_Std_Dev")
rnaseq_mapped_genes[,81]<-(rnaseq_mapped_genes[,75]-rnaseq_mapped_genes[,79])/rnaseq_mapped_genes[,80]
rnaseq_mapped_genes[,82]<-(rnaseq_mapped_genes[,76]-rnaseq_mapped_genes[,79])/rnaseq_mapped_genes[,80]
rnaseq_mapped_genes[,83]<-(rnaseq_mapped_genes[,77]-rnaseq_mapped_genes[,79])/rnaseq_mapped_genes[,80]
rnaseq_mapped_genes[,84]<-(rnaseq_mapped_genes[,78]-rnaseq_mapped_genes[,79])/rnaseq_mapped_genes[,80]
colnames(rnaseq_mapped_genes)[81:84]<-c("D0_zscore","D3_zscore","D6_zscore","D10_zscore")
# Fold changes
rnaseq_mapped_genes[,85]<-rnaseq_mapped_genes[,76]/rnaseq_mapped_genes[,75]#D3/D0
rnaseq_mapped_genes[,86]<-rnaseq_mapped_genes[,77]/rnaseq_mapped_genes[,75]#D6/D0
rnaseq_mapped_genes[,87]<-rnaseq_mapped_genes[,78]/rnaseq_mapped_genes[,75]#D10/D0
rnaseq_mapped_genes[,88]<-rnaseq_mapped_genes[,77]/rnaseq_mapped_genes[,76]#D6/D3
rnaseq_mapped_genes[,89]<-rnaseq_mapped_genes[,78]/rnaseq_mapped_genes[,76]#D10/D3
rnaseq_mapped_genes[,90]<-rnaseq_mapped_genes[,78]/rnaseq_mapped_genes[,77]#D10/D6
#read int the genelists
candidate_genes <- read.csv("~/Desktop/bibhas/candidate_genes.csv")
candidate_genes[,1] <- as.character(candidate_genes[,1])
Me<-subset(rnaseq_mapped_genes,rnaseq_mapped_genes[,1] %in% candidate_genes[,1])
candidate_genes[,2] <- as.character(candidate_genes[,2])
ES<-subset(rnaseq_mapped_genes,rnaseq_mapped_genes[,1] %in% candidate_genes[,2])
candidate_genes[,3] <- as.character(candidate_genes[,3])
iPSC<-subset(rnaseq_mapped_genes,rnaseq_mapped_genes[,1] %in% candidate_genes[,3])
```
First we see a decrease in the expression of Mesenchymal genes.
```{r message=FALSE, warning=FALSE}
library("gplots")
par(font.axis=2,font.lab=2,mfrow=c(1,3))
x<-heatmap.2(as.matrix(Me[,75:78]), col=greenred(75),
scale="row",key=T, symkey=T, Rowv=F,Colv=FALSE,cexRow=0.9,
trace='none',dendrogram="none",cexCol=1.2,srtCol =90,key.title=" ",key.xlab="",cex.axis=0.7,
ylab = "",labRow = Me[,1],labCol=c("3 Hour","3 Days","6 Days","10 Days"),density.info="none",
main="Mesenchymal genes")
me_zscore<-as.data.frame(x$carpet)
par(font.axis=2,font.lab=2,mfrow=c(1,2))
t<-rowSums(me_zscore)
t_up<-max(t)+5
t_dw<-min(t)-5
barplot(t,las=2,ylab="Summed zscore",names =c("3 Hour","3 Days","6 Days","10 Days"),ylim=c(t_dw,t_up))
box()
abline(h=0)
x<-colSums(Me[,75:78])
x_up<-max(x)+200
barplot(x,las=2,ylab="Summed Tpkm",names =c("3 Hour","3 Days","6 Days","10 Days"),ylim=c(0,x_up))
box()
```
Then we see an increase in the expression of Embryonic Stem cell (ES) genes.
```{r ES cells, message=FALSE, warning=FALSE}
par(font.axis=2,font.lab=2,mfrow=c(1,1))
x<-heatmap.2(as.matrix(ES[,75:78]), col=greenred(75),
scale="row",key=T, symkey=T, Rowv=F,Colv=FALSE,cexRow=0.9,
trace='none',dendrogram="none",cexCol=1.2,srtCol =90,key.title=" ",key.xlab="",cex.axis=0.7,
ylab = "",labRow = ES[,1],labCol=c("3 Hour","3 Days","6 Days","10 Days"),density.info="none",
main="ES genes")
es_zscore<-as.data.frame(x$carpet)
par(font.axis=2,font.lab=2,mfrow=c(1,2))
t<-rowSums(es_zscore)
t_up<-max(t)+5
t_dw<-min(t)-5
barplot(t,las=2,ylab="Summed zscore",names =c("3 Hour","3 Days","6 Days","10 Days"),ylim=c(t_dw,t_up))
box()
abline(h=0)
par(font.axis=2,font.lab=2)
x<-colSums(ES[,75:78])
x_up<-max(x)+200
barplot(x,las=2,ylab="Summed Tpkm",names =c("3 Hour","3 Days","6 Days","10 Days"),ylim=c(0,x_up))
box()
```
Also we see an increase in the expression of iPSC genes.
```{r iPSC cells, message=FALSE, warning=FALSE}
par(font.axis=2,font.lab=2,mfrow=c(1,1))
x<-heatmap.2(as.matrix(iPSC[,75:78]), col=greenred(75),
scale="row",key=T, symkey=T, Rowv=F,Colv=FALSE,cexRow=0.9,
trace='none',dendrogram="none",cexCol=1.2,srtCol =90,key.title=" ",key.xlab="",cex.axis=0.7,
ylab = "",labRow = iPSC[,1],labCol=c("3 Hour","3 Days","6 Days","10 Days"),density.info="none",
main="iPSC genes")
iPSC_zscore<-as.data.frame(x$carpet)
par(font.axis=2,font.lab=2,mfrow=c(1,2))
t<-rowSums(iPSC_zscore, na.rm = T)
t_up<-max(t)+5
t_dw<-min(t)-5
barplot(t,las=2,ylab="Summed zscore",names =c("3 Hour","3 Days","6 Days","10 Days"),ylim=c(t_dw,t_up))
box()
abline(h=0)
x<-colSums(iPSC[,75:78])
x_up<-max(x)+200
barplot(x,las=2,ylab="Summed Tpkm",names =c("3 Hour","3 Days","6 Days","10 Days"),ylim=c(0,x_up))
box()
```
Overall we find an increase in the stem cell marker genes and a decrease in the mesenchymal genes
```{r Fold change of candidate genese, fig.height=3, fig.width=3, message=FALSE, warning=FALSE}
# fold change
par(font.axis=2,font.lab=2)
Me[,87]<-Me[,78]/Me[,75]
ES[,87]<-ES[,78]/ES[,75]
iPSC[,87]<-iPSC[,78]/iPSC[,75]
boxplot(log(Me[,87],base=2),log(ES[,87],base=2),log(iPSC[,87],base=2),na.rm=T, las=1, col="gray",pch=18, names=c("Mesenchymal","ES","iPSC"),
ylab="log(2) Fold Change ( 10 Days/ 3hour)",lty=1)
abline(h=0, lwd=1,col="red")
```
We specifically see a monotonous decay in the expression of mesenchymal genes whereas we see a biphasic expression of stem cell genes
```{r Heatmap candidate genes,message=FALSE, warning=FALSE }
temp<-as.matrix(rbind(colMeans(Me[,75:78]),colMeans(ES[,75:78]),colMeans(iPSC[,75:78])))
colnames(temp)<-c("3 Hour","3 Days","6 Days","10 Days")
rownames(temp)<-c("Menchymal","ES","iPSC")
library("corrplot")
corrplot(as.matrix(temp),is.corr=F, tl.col="black",cl.pos="b", cl.ratio=0.3,method="color",addgrid.col="black")
```
Lastly, plot the bar plot of some important genes. Error bars are 95% CI.
```{r barplot_candidate_genes, fig.height=10, fig.width=8, message=FALSE, warning=FALSE}
imp<-c("Cdh2","Col1a1","Fn1","Snai1","Vim","Twist1","Col1a2","Alpl","Nanog","Pou5f1","Lin28a","Sox2","Fut4","Klf4","Lmna","Pcna")
imp_genes<-subset(rnaseq_mapped_genes,rnaseq_mapped_genes[,1] %in% imp)
error.bar <- function(x, y, upper, lower=upper, length=0.02,...){
if(length(x) != length(y) | length(y) !=length(lower) | length(lower) != length(upper))
stop("vectors must be same length")
arrows(x,y+upper, x, y-lower, angle=90, code=3, length=length, ...)
}
par(font.axis=2,font.lab=2,mfrow=c(4,4))
for( i in 1:nrow(imp_genes)){
y.means<-as.vector(imp_genes[i,c(75:78)])
a<-apply(imp_genes[i,c(63,67,71)],1,sd)
b<-apply(imp_genes[i,c(64,68,72)],1,sd)
c<-apply(imp_genes[i,c(65,69,73)],1,sd)
d<-apply(imp_genes[i,c(66,70,74)],1,sd)
y.sd<-(as.vector(c(a,b,c,d)))
ylimitup<-max(y.means)+1.5*max(y.sd)
barx <- barplot(as.matrix(y.means), names=c("3 Hour","3 Days","6 Days","10 Days"),main=imp_genes[i,1],ylim=c(0,ylimitup),las=2,ylab="TPkM",cex.axis = 1,
xpd = FALSE,cex.names=1)
error.bar(barx,as.matrix(y.means), (1.96*y.sd/sqrt(3)))
box()
}
```
#Session Information
```{r}
sessionInfo()
```