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//RtmpTk3QGD/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,] 37 20 78 57 62 4 37 124 37 37 37
## [2,] 66 87 47 75 100 31 66 200 66 66 66
## [3,] 65 6 85 32 142 48 65 284 65 65 65
## [4,] 48 18 106 97 84 9 48 168 48 48 48
## [5,] 38 51 41 84 112 7 38 224 38 38 38
## [6,] 91 30 89 73 85 15 91 170 91 91 91
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
## 12
a_j.match <-match(x[,"a"], x[,"j"], nomatch = 0)>0
table(a_j.match)["TRUE"]
## TRUE
## 43
(a_j.match_Per <- (table(a_j.match)["TRUE"]/nrow(x))*100)
## TRUE
## 86
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
## 2444 2851 3544 2915 5148 1310 2871 10296 3220 2698 3007
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] 1.883 1.495 1.163 1.046 0.964 ...
## $ rotation: num [1:11, 1:11] 0.125 0.572 0.147 0.132 -0.112 ...
## ..- 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] -1.872 1.088 -0.674 -1.06 -1.419 ...
## ..- 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,] -1.8724473 -0.10423090 -0.13109152 -1.7800161 -0.171980832
## [2,] 1.0881689 -0.76501224 0.14792732 0.1027551 0.013017462
## [3,] -0.6735578 0.34086219 1.22982930 1.9881347 1.165801153
## [4,] -1.0603338 0.34538218 0.68958982 -1.6172267 -0.105482707
## [5,] -1.4187249 -0.74352057 -0.37214221 -1.1970807 -1.078646538
## [6,] 1.4281097 -0.03368702 2.31011795 -1.0667025 0.556665611
## [7,] 0.6028013 -0.08226043 0.32048730 0.6723045 0.833314261
## [8,] -2.8911605 -0.09633714 -0.56728852 -0.0239698 -0.384765816
## [9,] 0.7784790 -0.14248240 -1.02853578 1.3843094 0.861718881
## [10,] 0.4637743 -0.54558844 -0.62304480 -0.9348074 0.006153661
## [11,] -0.3986697 2.67190582 -0.95025039 -1.4359945 1.160682926
## [12,] -2.3488969 1.42323525 -0.70316723 -0.8375198 1.465026283
## [13,] -0.2173296 1.04698869 -0.37728072 -0.2233423 0.842833946
## [14,] -2.6112801 -0.63733140 -0.55303274 0.1826531 -0.696609317
## [15,] -2.4744689 1.14959594 -0.64700158 0.7403774 0.957854278
## [16,] 2.1074577 -2.80109083 1.38828815 0.2080810 -0.063636392
## [17,] 1.9719846 -0.76510462 0.65744158 -1.3173462 -1.261644881
## [18,] 0.8170270 -0.57174218 0.42702213 0.8119207 1.483901305
## [19,] 1.7139877 0.20025177 0.92138819 -0.2204587 0.900005840
## [20,] -2.1005643 -3.05102788 0.98803206 1.1772021 0.093879492
## [21,] -0.9438407 -2.29491431 -0.81818699 0.7224788 0.417841141
## [22,] -2.3973854 -2.22482910 -1.27477136 -0.3410124 -0.107845111
## [23,] 0.2365796 -3.64453693 0.10016992 -1.9442957 -0.245682020
## [24,] 0.1407406 2.84307502 -1.19064463 -0.2639410 -0.579441503
## [25,] -4.1364119 -0.10189934 -1.29216113 1.1781582 -0.631014571
## [26,] -0.4970310 -0.73335898 0.78900322 1.3382432 1.061190601
## [27,] -0.7795006 -1.18348022 -0.10888800 -0.9593339 -0.297396954
## [28,] 2.2073764 1.84371005 -0.96821145 0.3488827 0.501314483
## [29,] 2.8946159 1.78972432 -0.58957689 -1.4896832 -0.172162585
## [30,] 0.7360835 -2.18824687 -1.26500635 1.3965849 -0.163021994
## [31,] -1.2425495 0.11241876 -0.65272526 -1.5660168 0.162229409
## [32,] -1.5183155 2.34580948 -0.01127993 0.0544759 1.397377967
## [33,] 0.6074943 -0.08962740 -0.26187598 0.2799208 2.561802005
## [34,] 3.6575264 -0.76870863 0.99462060 -0.9173675 1.419102536
## [35,] -0.9732570 0.18226939 2.67354773 -0.7118090 -0.847108471
## [36,] 2.5044734 1.53128954 -2.76190181 -0.8286679 -1.219736392
## [37,] 0.8292500 -2.31603834 -2.38636742 -0.7537585 -1.055061413
## [38,] -0.3412793 -1.85589475 0.88553920 -0.7190147 0.509838874
## [39,] -2.6690295 -0.20742961 -0.65317899 0.7972965 0.452073332
## [40,] -1.5264470 1.81083535 2.25196658 -0.9196952 0.015329800
## [41,] -1.3012529 1.94916465 0.75233248 0.5652313 -1.742334253
## [42,] 1.6352391 0.69860462 -1.01973282 0.2041362 -0.098840681
## [43,] 3.5296939 0.21355277 -1.27375429 1.3707425 -0.304467155
## [44,] -1.0230182 2.03998413 1.20822031 0.8872705 -0.402859760
## [45,] 2.3252442 -0.17910134 0.49066698 1.7546107 -1.249389662
## [46,] 0.6384733 0.71975676 2.17650605 0.8227703 -1.581948810
## [47,] 1.1751934 1.08651655 1.67194859 -0.1346704 -0.850328760
## [48,] 3.6937299 -0.14435981 -0.62501258 0.5465019 -0.356070562
## [49,] -1.9917060 1.06856618 0.35095076 1.0792217 -2.191164021
## [50,] 1.6249544 0.85834223 -0.31948487 1.5894669 -0.980314084
## PC6 PC7 PC8 PC9 PC10
## [1,] 0.05039552 -0.008870074 1.08020103 -0.568931825 -0.30775785
## [2,] -0.66697429 0.184434512 -0.46257354 0.533722400 0.26231667
## [3,] 0.59695280 -0.319490855 -0.36012379 -0.893957967 -0.10951585
## [4,] -0.41546225 -1.133911128 0.59366986 -0.311891739 -0.96920051
## [5,] -0.38824239 -0.528684873 -0.20423108 0.433643359 0.04379098
## [6,] 0.25024797 -0.493228051 -0.10034949 -1.079131393 0.20343191
## [7,] 0.67102964 0.465994924 -0.17871298 0.038001000 -0.11496641
## [8,] 0.68802194 -0.365217369 0.06545451 -0.050346668 -0.33077769
## [9,] -0.23194571 0.704470641 -0.11185925 1.211739711 -0.68587950
## [10,] 1.31511029 0.975862648 -0.08813331 0.667987246 0.05648386
## [11,] -0.20996857 0.446271622 0.24936280 0.202751126 0.38195601
## [12,] 0.04732737 -0.381729323 -0.07112440 -0.385168470 0.41278502
## [13,] -1.90609712 -0.181426899 0.42532329 0.152663159 -0.11763646
## [14,] -0.19693309 -0.531004858 -0.37208191 0.128316813 0.19109232
## [15,] -0.23192452 -0.977209524 -0.86130141 0.166432746 0.13290418
## [16,] -0.90141894 0.706668111 0.65843808 -0.419689455 -0.09504969
## [17,] 0.51814942 0.586522031 -0.22425654 0.566342470 0.16336929
## [18,] -0.64015041 0.714999848 0.46480485 -0.470168930 0.10614930
## [19,] -0.27968529 -0.311183276 -0.54216342 0.199521743 -0.07631427
## [20,] -0.54859434 -1.344407031 -0.57825723 -0.642159955 -0.10511415
## [21,] -0.73015376 0.663389641 0.35937209 0.182240768 0.08681074
## [22,] -0.04408387 1.062695045 0.83855003 -0.237719833 0.47360379
## [23,] -0.77097909 -0.069555914 0.35111981 0.209525349 -0.06605614
## [24,] -1.35420595 0.554919645 0.23146709 1.027491049 0.05092978
## [25,] -0.45743353 0.043846543 0.19969762 -0.220684253 0.37966925
## [26,] 0.13707838 0.729001424 0.14846260 0.147303811 0.40064359
## [27,] -0.39902180 0.096898568 0.28519843 -0.062300315 0.17947712
## [28,] -0.14805386 0.881745579 -0.07985853 -0.765616006 -0.04948504
## [29,] 0.26868651 0.510510156 -0.24335288 -0.731128184 -0.04252973
## [30,] -1.05796018 -0.715356712 -1.18914162 0.348970854 -0.64752565
## [31,] 1.40869085 -0.033928200 -0.59230307 -0.082627128 0.49408581
## [32,] 0.74430851 -0.604622571 -0.96142202 0.301994406 0.41376292
## [33,] -0.02413439 0.742970373 0.12365185 -0.143053250 0.23286596
## [34,] -0.24384151 0.223932975 -0.29193708 0.096860384 -0.18153399
## [35,] -0.69674875 -0.343046167 0.52957697 0.498779868 0.43259227
## [36,] -1.16411490 -0.285661291 -0.86120735 -0.893384461 -0.54889682
## [37,] 1.40727145 0.243721577 -1.30617702 -0.385850057 -0.07645911
## [38,] 1.75027146 0.076244897 -0.15868050 0.515040823 -0.32260942
## [39,] 0.90575176 0.159432047 0.37533884 0.156510778 -0.59089522
## [40,] 0.40984575 -0.900372530 0.16553755 0.719959959 -0.20309509
## [41,] -0.70891800 0.310803931 -0.22868832 -0.268162881 0.46336090
## [42,] 0.57641262 -0.916225762 1.36202763 0.004323762 -0.16698472
## [43,] 1.24477581 -1.397412436 1.44606256 0.207810924 -0.06369769
## [44,] -0.20406680 0.807292594 0.13249437 -0.097480774 -1.09909972
## [45,] 0.13170586 -0.118020367 -0.14513688 0.047332800 0.75046227
## [46,] 0.14610126 0.810480109 -1.00811691 0.040597280 -0.04364749
## [47,] 0.06697538 0.019844331 -0.65239816 0.177953404 -0.30712422
## [48,] -0.33998529 -1.583152080 0.29771336 0.225224239 1.12315852
## [49,] 0.49257574 1.198121906 0.52510612 -0.528533831 0.22612598
## [50,] 1.13341232 -0.377358384 0.96495736 0.028945146 -0.33997600
## PC11
## [1,] -1.837111e-17
## [2,] -1.635384e-16
## [3,] 1.362507e-16
## [4,] 4.309512e-16
## [5,] -1.099934e-16
## [6,] 2.798281e-16
## [7,] -2.476053e-17
## [8,] -7.911540e-17
## [9,] -1.492270e-16
## [10,] -1.632493e-16
## [11,] 1.648590e-16
## [12,] -1.265540e-16
## [13,] 4.923436e-17
## [14,] -2.528826e-16
## [15,] -4.795647e-17
## [16,] -1.361398e-16
## [17,] 9.169083e-17
## [18,] -1.051478e-16
## [19,] 1.162477e-16
## [20,] -3.324023e-16
## [21,] -5.232438e-16
## [22,] -5.997701e-16
## [23,] -3.985177e-16
## [24,] 2.837481e-16
## [25,] -3.824239e-16
## [26,] -1.799801e-16
## [27,] -2.453175e-16
## [28,] 4.659766e-17
## [29,] 2.037117e-16
## [30,] -5.581404e-16
## [31,] -1.963944e-16
## [32,] 1.332205e-16
## [33,] -2.157811e-16
## [34,] -6.054974e-19
## [35,] 4.765937e-16
## [36,] -1.371786e-16
## [37,] -6.993774e-16
## [38,] -5.921315e-17
## [39,] -1.281948e-16
## [40,] 7.194131e-16
## [41,] 3.153056e-16
## [42,] 3.557245e-16
## [43,] 4.006007e-16
## [44,] 4.999865e-16
## [45,] 9.271241e-17
## [46,] 2.537424e-16
## [47,] 4.210089e-16
## [48,] 9.591185e-17
## [49,] 1.991511e-16
## [50,] 4.612755e-16
sum(prcomp$sdev^2)/
prcomp$sdev[1]^2
## [1] 3.096165
(PC1.Perc <- (prcomp$sdev[1]^2/sum(prcomp$sdev^2))*100)
## [1] 32.29802
(PC2.Perc <- (prcomp$sdev[2]^2/sum(prcomp$sdev^2))*100)
## [1] 20.35297
(PC3.Perc <- (prcomp$sdev[3]^2/sum(prcomp$sdev^2))*100)
## [1] 12.31197
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] 1.864 1.48 1.151 1.036 0.954 ...
## ..- attr(*, "names")= chr [1:11] "Comp.1" "Comp.2" "Comp.3" "Comp.4" ...
## $ loadings: 'loadings' num [1:11, 1:11] 0.125 0.572 0.147 0.132 -0.112 ...
## ..- 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] -1.872 1.088 -0.674 -1.06 -1.419 ...
## ..- 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
## 1.8640861 1.4797621 1.1509113 1.0356927 0.9539285 0.7574435 0.6776082
## Comp.8 Comp.9 Comp.10 Comp.11
## 0.6021938 0.4760694 0.4058410 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.125 0.302 0.672 0.116 0.258 0.398 0.365 0.251
## b 0.572 0.159 -0.461 0.271 0.141 0.568
## c 0.147 -0.501 0.133 0.335 -0.329 -0.186 -0.190 0.115 0.636
## d 0.132 0.161 0.200 -0.238 0.685 -0.546 0.185 0.220
## e -0.112 -0.271 -0.171 -0.302 -0.267 0.430 0.179
## f 0.108 -0.822 0.365 0.387
## g 0.507 0.279 -0.190 -0.119 -0.240 0.210 -0.696 0.165
## h -0.112 -0.271 -0.171 -0.302 -0.267 0.430 0.179
## i 0.389 -0.421 0.359 -0.209 -0.622 0.255 -0.158 -0.134
## j 0.343 0.295 0.361 -0.101 0.200 -0.247 -0.402 -0.442 0.238 -0.374
## k 0.235 -0.530 0.212 0.371 0.406 -0.550
## 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
## 32.29802
(Comp2.Perc <- (princomp$sdev[2]^2/sum(princomp$sdev^2))*100)
## Comp.2
## 20.35297
(Comp3.Perc <- (princomp$sdev[3]^2/sum(princomp$sdev^2))*100)
## Comp.3
## 12.31197
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")
