if (!requireNamespace("BiocManager"))
install.packages("BiocManager")
## Loading required namespace: BiocManager
BiocManager::install("gplots", dependencie=T)
## Bioconductor version 3.8 (BiocManager 1.30.4), R 3.5.2 (2018-12-20)
## Installing package(s) 'gplots'
##
## The downloaded binary packages are in
## /var/folders/1s/2pw7wx312qg5vcs25c1qf3p00000gn/T//RtmpDejZlQ/downloaded_packages
library(gplots)
##
## Attaching package: 'gplots'
## The following object is masked from 'package:stats':
##
## lowess
a<- sample(1:100, 50, replace = T)
b<- sample(1:120, 50, replace=T)
c<- sample(1:125, 50, replace=F)
d<- sample(20:100, 50, replace=F)
e<- sample(50:150, 50, replace= T)
f<- sample(1:50, 50, replace=T)
g<- c(a[1:25], b[26:50]) # g has half value matching a and half values matching b
h<- c(2*e) # h has values similar to e but twice of its e value
i<- c(a[1:10], b[11:20],c[21:30], d[31:40], e[41:50]) # values of i are frpm a to e samples
j<- c(a[1:40],b[41:50]) # j values are 80% matching with a and 20% matching with b
k<- c(a[1:10], c[11:35], d[36:50]) # k values 20%, 50% and 30% matchng with a,c and d respectively
x<- cbind(a,b,c,d,e,f,g,h,i,j,k)
head(x)
## a b c d e f g h i j k
## [1,] 29 75 10 41 77 5 29 154 29 29 29
## [2,] 6 32 15 55 141 8 6 282 6 6 6
## [3,] 82 20 69 47 140 13 82 280 82 82 82
## [4,] 69 21 50 71 85 38 69 170 69 69 69
## [5,] 87 59 74 59 119 16 87 238 87 87 87
## [6,] 13 91 94 84 60 3 13 120 13 13 13
colnames(x)
## [1] "a" "b" "c" "d" "e" "f" "g" "h" "i" "j" "k"
a_b.match <- match(x[,"a"], x[,"b"], nomatch = 0)>0
table(a_b.match)["TRUE"]
## TRUE
## 14
a_j.match <-match(x[,"a"], x[,"j"], nomatch = 0)>0
table(a_j.match)["TRUE"]
## TRUE
## 44
(a_j.match_Per <- (table(a_j.match)["TRUE"]/nrow(x))*100)
## TRUE
## 88
rc<- rainbow(nrow(x), start = 0, end= 0.3)
cc <- rainbow(ncol(x), start = 0, end = .3)
hv <- heatmap(x, col = cm.colors(256), scale = "column",
RowSideColors = rc, ColSideColors = cc, margins = c(5,10),
xlab = "Sample", ylab = "Gene",
main= "raw values")
colfunc<-colorRampPalette(c("red","yellow","springgreen"))
hv <- heatmap(x, col=c(colfunc(ncol(x))), scale = "column",
RowSideColors = c(colfunc(nrow(x))), ColSideColors = c(colfunc(ncol(x))),
margins = c(5,10),
xlab = "Sample", ylab = "Gene",
main= "raw values with color grad")
colSums(x)
## a b c d e f g h i j k
## 2354 3264 3141 2860 4960 1187 2699 9920 3685 2500 3149
x1<- colSums(x)
x2<- cbind(a=c(x[,1]/x1[1]), b=c(x[,2]/x1[2]),c=c(x[,3]/x1[3]),d=c(x[,4]/x1[4]),
e=c(x[,5]/x1[5]),f=c(x[,6]/x1[6]),g=c(x[,7]/x1[7]),h=c(x[,8]/x1[8]),
i=c(x[,9]/x1[9]), j=c(x[,10]/x1[10]), k=c(x[,11]/x1[11]))
# Optional
#write.csv(x, file= "Array card_Raw.csv")
#write.csv(x2, file= "Array card Normalized.csv")
#write.table(x2,file= "mydata.txt", sep="\t")
colfunc<-colorRampPalette(c("springgreen","yellow","red"))
hv <- heatmap(x2, col=c(colfunc(ncol(x2))), scale = "column",
RowSideColors = c(colfunc(nrow(x2))), ColSideColors = c(colfunc(ncol(x2))),
xlab = "Sample", ylab = "Gene",
main= "Normalized values"
)
#for color scale in heat map
heatmap.2(x2, col=c(colfunc(ncol(x2))), scale = "column",
RowSideColors = c(colfunc(nrow(x2))), ColSideColors = c(colfunc(ncol(x2))),
xlab = "Sample", ylab = "Gene",
main= "Normalized values_color scale",
density.info = "none",
trace= "none")
# Normalized per 100 unit
x3<- x2*100
heatmap.2(x3, col=c(colfunc(ncol(x3))), scale = "column",
RowSideColors = c(colfunc(nrow(x3))), ColSideColors = c(colfunc(ncol(x3))),
xlab = "Sample", ylab = "Gene",
main= "Normalized values per 100_color scale",
density.info = "none",
trace= "none")
prcomp<-prcomp(x3)
str(prcomp)
## List of 5
## $ sdev : num [1:11] 2.12 1.55 1.09 1.05 0.95 ...
## $ rotation: num [1:11, 1:11] 0.478 0.209 0.223 0.101 0.025 ...
## ..- attr(*, "dimnames")=List of 2
## .. ..$ : chr [1:11] "a" "b" "c" "d" ...
## .. ..$ : chr [1:11] "PC1" "PC2" "PC3" "PC4" ...
## $ center : Named num [1:11] 2 2 2 2 2 2 2 2 2 2 ...
## ..- attr(*, "names")= chr [1:11] "a" "b" "c" "d" ...
## $ scale : logi FALSE
## $ x : num [1:50, 1:11] -2.322 -4.176 1.85 0.838 2.542 ...
## ..- attr(*, "dimnames")=List of 2
## .. ..$ : NULL
## .. ..$ : chr [1:11] "PC1" "PC2" "PC3" "PC4" ...
## - attr(*, "class")= chr "prcomp"
prcomp$x
## PC1 PC2 PC3 PC4 PC5
## [1,] -2.32249286 0.32887916 -1.111437020 1.29770628 1.15380228
## [2,] -4.17601566 0.39018043 0.065156075 1.23792662 1.29957751
## [3,] 1.85014902 -0.34026317 1.081256188 -0.05649963 2.04423394
## [4,] 0.83784485 1.04537546 0.923578592 -1.32189829 0.61976054
## [5,] 2.54214565 -0.19420439 0.433505568 0.13399596 1.14717151
## [6,] -2.70228619 -1.33106697 -2.098313558 0.67592540 -0.23130729
## [7,] -3.56906336 2.14560894 -0.918570604 0.37470802 -1.95419479
## [8,] -2.67882910 1.42367550 -0.412159832 -0.39059807 -0.73568049
## [9,] 3.11828238 -1.12744331 -0.447397298 -1.44648584 0.93175743
## [10,] 3.73623764 -0.06499904 0.252339914 -0.07049437 -0.56858318
## [11,] 0.27831527 -2.81850845 0.128370276 -0.18183730 0.51011038
## [12,] -1.86251752 0.27262186 -0.548117839 -0.28553896 -0.63307997
## [13,] 4.17900647 -0.44824052 0.265879194 -0.18713518 -1.15342189
## [14,] -2.76965219 2.48881028 -0.245711351 -1.79564495 1.47138248
## [15,] 0.11882056 -1.61473488 -0.726877848 0.59255245 -0.59918867
## [16,] -1.13094973 0.67637553 0.279508424 1.42203807 -0.71837485
## [17,] -2.78627771 -1.42198667 -0.334225102 1.93351056 -0.77729247
## [18,] -2.14079134 -1.85339896 0.574065644 -0.28691990 0.05578737
## [19,] 0.88289247 -0.86330992 -0.629558915 0.88685326 0.52380862
## [20,] -1.29731494 -0.72083585 -0.173027181 0.30818935 -1.47531949
## [21,] 0.67659759 0.97881666 -0.631786844 -0.90340205 -2.43371905
## [22,] 2.46205594 0.87523082 -0.431340303 0.75036722 -0.33789483
## [23,] 1.93531296 -0.29298513 1.067106351 0.03538100 0.59448990
## [24,] -1.42736376 -0.00532006 -0.008501324 0.28263667 1.34537587
## [25,] 1.26588133 -1.91883796 0.299972059 -0.74217625 -0.09295698
## [26,] 3.09835485 -1.56174697 -0.053208537 0.60519089 -0.54078416
## [27,] -0.07961498 -2.36903229 -0.603060899 -0.73428154 -1.47890377
## [28,] -1.07258286 -1.16943492 -1.111394121 0.72992538 0.59764744
## [29,] 4.15101848 0.16495074 -1.036842097 -0.90704473 -0.07681942
## [30,] 1.14496491 -2.72826065 0.951685148 -0.28259812 -0.60100216
## [31,] -0.76375505 -1.53708707 -0.785040563 -0.83193357 0.86279950
## [32,] 2.31215375 3.94828739 -1.434889690 -0.29590151 0.59637551
## [33,] 0.18997718 3.27293780 0.262938017 -0.03552056 0.47430765
## [34,] -2.53604433 -2.46456074 0.477789304 -1.26031937 0.02473195
## [35,] 1.08014360 2.97792462 -0.389146621 0.26803913 -0.78080090
## [36,] 0.21491264 1.06484254 1.260148528 0.75529193 0.21781175
## [37,] -0.01939498 0.16119828 -1.525129741 -2.66875485 0.25705575
## [38,] 2.50562777 1.97014431 0.912412619 0.21377513 -1.21204505
## [39,] -1.91916610 -0.88025204 1.553184323 -0.79827987 0.04034540
## [40,] 0.06349961 2.17767808 0.615604327 -1.23530167 -0.27959014
## [41,] -1.38059887 1.36158120 1.126024957 0.23638707 -1.07141911
## [42,] 1.01156865 -0.72198029 1.451966113 0.92815658 0.47293060
## [43,] 1.20361312 0.47665410 -2.243294904 0.71836931 1.41243325
## [44,] -1.16914903 0.21439954 3.837087846 -0.35304058 -0.61477114
## [45,] -2.20490739 -1.31958878 -0.827544269 -2.21365312 -0.58478237
## [46,] -1.82731156 0.53947575 0.722447474 -0.83200734 1.51473596
## [47,] -2.20460631 1.17861422 0.341427255 0.07625829 0.21936433
## [48,] -0.30625829 0.15931317 1.503289302 1.41687759 0.30969649
## [49,] 2.08500256 -0.64649156 -2.149356640 1.33697852 0.04353556
## [50,] 1.40256488 0.12099421 0.489189603 2.90022697 0.21090320
## PC6 PC7 PC8 PC9 PC10
## [1,] -0.20242196 0.66187009 -0.416729702 -0.679508397 0.226602110
## [2,] 1.21900153 -0.26133240 0.224928989 0.002303584 -0.270067694
## [3,] 0.46806079 -0.64373350 -0.137821480 0.202706640 -0.105073926
## [4,] -0.36121402 0.15360976 0.781798567 -0.330651688 0.016398405
## [5,] 0.36910159 -0.04368570 -0.146213009 -0.033030391 0.062948845
## [6,] 1.05739959 1.55507915 0.993554173 0.241114053 -0.204207497
## [7,] 1.82994243 -0.71222750 -0.724761821 0.098760397 0.063970726
## [8,] 1.25463444 -0.89890093 -0.549984295 0.284244035 -0.026479032
## [9,] -0.45240413 0.58520699 0.312446184 0.217312862 0.072420816
## [10,] 1.19348178 0.38912726 0.433623619 0.170194801 0.002010667
## [11,] 1.11932021 -0.40428036 -0.177040159 -0.588311725 0.196102603
## [12,] -0.97168499 0.53661355 -0.192417953 -0.065826934 0.241729036
## [13,] 0.80141415 0.04196627 -0.306960225 0.122797822 0.236783687
## [14,] -0.46537973 -0.40590042 0.176994266 -0.335940008 -0.024225750
## [15,] -0.55192531 0.23069211 -1.119056748 0.122176168 0.396105432
## [16,] 0.22962780 0.32898227 -0.124213742 0.538475708 -0.042833052
## [17,] -0.38970616 0.50936095 -0.668425619 0.361840360 0.167358403
## [18,] 1.14612476 -0.43362699 0.373925148 -0.563319700 0.039875026
## [19,] -0.16530908 0.28987129 -0.826792168 0.159514082 0.217938106
## [20,] -0.53334327 1.20489542 0.591944716 -0.215781443 0.195573816
## [21,] -0.10926024 0.89163336 0.293447890 -0.500949507 0.384436283
## [22,] 1.24578027 0.38551046 -0.128830702 -0.153099494 0.119958788
## [23,] 0.91069771 0.15434230 0.809961225 0.109356144 -0.168617752
## [24,] 0.18340969 0.29671507 0.386831694 -0.192004323 -0.076159470
## [25,] -0.77146003 0.65729489 0.386092457 -0.095234202 0.195732391
## [26,] 0.24964973 0.21795905 0.170663201 -0.199142083 -0.503092018
## [27,] -1.58238769 -1.03220833 -0.742728007 -0.394994216 -1.012898821
## [28,] -0.79656936 0.34136692 0.188762189 -0.677973274 -0.807979813
## [29,] -0.35161380 -0.38109900 -1.332730957 -0.067751245 0.181543563
## [30,] -0.12364986 -0.50855549 -0.331364580 0.307481098 -0.024265506
## [31,] 0.33471998 -0.47607906 -0.446747712 -0.874470290 0.092192409
## [32,] 0.25632064 -0.11341226 -0.649212460 -0.206174417 -0.223530736
## [33,] 0.29276154 -0.01370216 0.011289772 0.502940000 -0.032337210
## [34,] 0.42988096 -0.92229928 -0.017683007 -0.633520609 0.039137691
## [35,] -0.61672870 0.52293863 0.641406826 0.056833861 -0.868763993
## [36,] 0.48711346 -0.35496358 -0.017893331 -0.537968779 0.019270115
## [37,] 0.57700616 -0.01374845 0.083435522 0.978888312 0.094415315
## [38,] 0.46476933 -0.42528540 0.120564330 -0.780467969 -0.367417158
## [39,] 0.51169309 -0.18086190 0.586965676 0.420422741 0.333644333
## [40,] -0.93496722 0.50732946 0.269124168 -0.873157074 0.619570341
## [41,] -1.43978203 -1.20907567 0.714605993 -0.202387025 0.322224332
## [42,] 0.04771048 0.70910282 -0.153072415 0.489771536 -0.120308047
## [43,] -0.80178369 -0.14133847 -0.002811591 0.339530779 0.198758692
## [44,] -0.20509838 0.51720036 -0.783736886 1.001399642 -0.320081991
## [45,] -0.39861458 -0.25160866 0.930055598 0.791061278 -0.120503475
## [46,] -0.69038368 1.42969622 -1.107009850 0.291264126 -0.143685003
## [47,] -1.52200379 -0.83093514 -0.208824599 0.552838328 0.108356388
## [48,] -0.80259601 0.01181899 0.105882863 -0.058655158 0.119370762
## [49,] -0.51149856 -1.60490957 0.950573511 1.150559477 0.146379683
## [50,] -0.92783589 -0.86641341 0.774184442 -0.253467888 0.351719179
## PC11
## [1,] -2.814181e-16
## [2,] 3.644871e-16
## [3,] 1.899766e-16
## [4,] -5.918947e-17
## [5,] 6.118945e-18
## [6,] -3.724267e-17
## [7,] 3.048034e-16
## [8,] 2.528307e-16
## [9,] -3.997833e-16
## [10,] 1.595007e-16
## [11,] -9.171420e-17
## [12,] -2.378483e-16
## [13,] -5.186127e-17
## [14,] -3.479819e-16
## [15,] -3.795839e-16
## [16,] 2.963897e-16
## [17,] 5.920528e-18
## [18,] 1.154093e-16
## [19,] -2.549532e-16
## [20,] -2.428613e-16
## [21,] -3.901106e-16
## [22,] 8.299056e-17
## [23,] 3.301299e-16
## [24,] 3.098402e-17
## [25,] -2.901988e-16
## [26,] 5.576773e-16
## [27,] 6.528988e-16
## [28,] 4.992605e-16
## [29,] -4.386043e-16
## [30,] 9.719966e-17
## [31,] -3.067886e-16
## [32,] -3.248997e-18
## [33,] 8.021044e-17
## [34,] -7.254190e-17
## [35,] 7.651838e-16
## [36,] 4.033121e-16
## [37,] -6.807279e-16
## [38,] 8.097148e-16
## [39,] -1.148993e-16
## [40,] -6.348460e-16
## [41,] -8.592281e-17
## [42,] 3.832430e-16
## [43,] -5.880533e-16
## [44,] 7.395310e-16
## [45,] -4.344822e-16
## [46,] -6.621912e-17
## [47,] -2.082638e-16
## [48,] 2.405153e-16
## [49,] -3.579730e-16
## [50,] 5.679174e-17
sum(prcomp$sdev^2)/
prcomp$sdev[1]^2
## [1] 2.630237
(PC1.Perc <- (prcomp$sdev[1]^2/sum(prcomp$sdev^2))*100)
## [1] 38.01939
(PC2.Perc <- (prcomp$sdev[2]^2/sum(prcomp$sdev^2))*100)
## [1] 20.36982
(PC3.Perc <- (prcomp$sdev[3]^2/sum(prcomp$sdev^2))*100)
## [1] 10.09124
plot(prcomp)
plot(prcomp$rotation,col= c(1:nrow(prcomp$rotation)), cex= 4, main= "PCA by prcomp",
xlab= paste("PC1", "[", round(PC1.Perc, digits = 1),"%", "]"),
ylab= paste("PC2", "[", round(PC2.Perc, digits = 1),"%", "]"))
text(prcomp$rotation, labels=row.names(prcomp$rotation))
bi_prcomp<- biplot(prcomp, main= "Bi-plot_prcomp")
princomp<- princomp(x3)
plot(princomp)
str(princomp)
## List of 7
## $ sdev : Named num [1:11] 2.1 1.54 1.08 1.04 0.94 ...
## ..- attr(*, "names")= chr [1:11] "Comp.1" "Comp.2" "Comp.3" "Comp.4" ...
## $ loadings: 'loadings' num [1:11, 1:11] 0.478 0.209 0.223 0.101 0.025 ...
## ..- attr(*, "dimnames")=List of 2
## .. ..$ : chr [1:11] "a" "b" "c" "d" ...
## .. ..$ : chr [1:11] "Comp.1" "Comp.2" "Comp.3" "Comp.4" ...
## $ center : Named num [1:11] 2 2 2 2 2 2 2 2 2 2 ...
## ..- attr(*, "names")= chr [1:11] "a" "b" "c" "d" ...
## $ scale : Named num [1:11] 1 1 1 1 1 1 1 1 1 1 ...
## ..- attr(*, "names")= chr [1:11] "a" "b" "c" "d" ...
## $ n.obs : int 50
## $ scores : num [1:50, 1:11] -2.322 -4.176 1.85 0.838 2.542 ...
## ..- attr(*, "dimnames")=List of 2
## .. ..$ : NULL
## .. ..$ : chr [1:11] "Comp.1" "Comp.2" "Comp.3" "Comp.4" ...
## $ call : language princomp(x = x3)
## - attr(*, "class")= chr "princomp"
princomp$sdev
## Comp.1 Comp.2 Comp.3 Comp.4 Comp.5 Comp.6 Comp.7
## 2.1028409 1.5392085 1.0833680 1.0424898 0.9399580 0.7932800 0.6541359
## Comp.8 Comp.9 Comp.10 Comp.11
## 0.5636795 0.4753038 0.3079254 0.0000000
princomp$n.obs
## [1] 50
princomp$loadings
##
## Loadings:
## Comp.1 Comp.2 Comp.3 Comp.4 Comp.5 Comp.6 Comp.7 Comp.8 Comp.9 Comp.10
## a 0.478 0.184 0.196 0.235 0.301 0.236 0.561 0.420
## b 0.209 0.114 -0.507 -0.425 -0.569 0.165 0.145 0.335 -0.120
## c 0.223 -0.603 -0.253 0.425 -0.210 0.317 -0.124 0.409 0.113
## d 0.101 0.291 -0.228 -0.273 0.258 0.448 -0.700 0.101
## e 0.293 -0.224 0.438 -0.352 0.126 0.147
## f 0.552 0.158 0.591 -0.530 -0.192
## g 0.514 0.135 -0.214 0.121 -0.181 -0.311 -0.168 -0.207 0.672
## h 0.293 -0.224 0.438 -0.352 0.126 0.147
## i 0.253 -0.139 0.404 -0.258 -0.265 -0.573 0.133 0.517
## j 0.509 0.175 -0.145 0.250 -0.224 -0.326 0.145 -0.668
## k 0.281 -0.452 0.360 -0.198 -0.163 0.147 -0.663 -0.222
## Comp.11
## a
## b
## c
## d
## e 0.707
## f
## g
## h -0.707
## i
## j
## k
##
## Comp.1 Comp.2 Comp.3 Comp.4 Comp.5 Comp.6 Comp.7 Comp.8
## SS loadings 1.000 1.000 1.000 1.000 1.000 1.000 1.000 1.000
## Proportion Var 0.091 0.091 0.091 0.091 0.091 0.091 0.091 0.091
## Cumulative Var 0.091 0.182 0.273 0.364 0.455 0.545 0.636 0.727
## Comp.9 Comp.10 Comp.11
## SS loadings 1.000 1.000 1.000
## Proportion Var 0.091 0.091 0.091
## Cumulative Var 0.818 0.909 1.000
# PCA_princomp contribution in %
(Comp1.Perc <- (princomp$sdev[1]^2/sum(princomp$sdev^2))*100)
## Comp.1
## 38.01939
(Comp2.Perc <- (princomp$sdev[2]^2/sum(princomp$sdev^2))*100)
## Comp.2
## 20.36982
(Comp3.Perc <- (princomp$sdev[3]^2/sum(princomp$sdev^2))*100)
## Comp.3
## 10.09124
plot(princomp$loadings, col= c(1:8), cex= 4, main= "PCA by Princomp",
xlab= paste("PC1", "[", round(Comp1.Perc, digits = 1),"%", "]"),
ylab= paste("PC2", "[", round(Comp2.Perc, digits = 1),"%", "]"))
text(princomp$loadings, labels=row.names(princomp$loadings))
biplot(princomp, main= "Bi-Plot PCA by princomp")
