• PCA


    #加载R包

    > library(psych)
    > library(reshape2)
    > library(ggplot2)
    > library(factoextra)

    ##提前把名为ehbio_salmon.DESeq2.normalized.symbol.txt的文件放在工作目录,以制表符分割的.txt文件

    > exprData <- "ehbio_salmon.DESeq2.normalized.symbol.txt"  ##这里只是定义exprData和sampleFile两个变量,文本类型的(一开始我还以为就把文件内容给赋值
    > sampleFile <- "sampleFile"  ##了)

    > data <- read.table(exprData,header = T,row.names = NULL,sep =" ")  ###read.table()读入文件。
    > class(data$id)  ##判断某一列的数据类型
    [1] "character

    > data <- read.table(exprData,header = T,row.names = "id",sep =" ")  ##一开始不明白row.names=NULL是怎么回事,并且我试了不加这个参数时data是一样的。
    Error in read.table(exprData, header = T, row.names = "id", sep = " ") :   ##看文档知道row.names是指定哪一列为行名的,
    'row.names'里不能有重复的名字   ##看了错误提示有些懵懂,之所以不能以id列为行名,是因为这列id有相同的。我下面验证了一下

    > rownames(data) <- data$id
    Error in `.rowNamesDF<-`(x, value = value) : 不允许有重复的'row.names'
    此外: Warning message:
    non-unique values when setting 'row.names': ‘ALG1L9P’, ‘ATXN7’, ‘BMS1P21’, ‘BMS1P4’, ‘CCDC39’, ‘CYB561D2’, ‘DIABLO’, ‘DNAJC9-AS1’, ‘DUXAP8’, ‘GOLGA8M’, ‘HSPA14’, ‘IPO5P1’, ‘ITFG2-AS1’, ‘LINC-PINT’, ‘LINC00484’, ‘LINC00941’, ‘LINC01238’, ‘LINC01297’, ‘LINC01422’, ‘LINC01481’, ‘MATR3’, ‘OR7E47P’, ‘PKD1P1’, ‘POLR2J3’, ‘POLR2J4’, ‘RAET1E-AS1’, ‘RF00012’, ‘RF00017’, ‘RF00019’, ‘RF02271’, ‘RGS5’, ‘RMRP’, ‘SCO2’, ‘SNHG28’, ‘SNX29P2’, ‘SPATA13’, ‘TBCE’, ‘TMSB15B’, ‘ZNF503’

    # 处理重复名字,谨慎处理,先找到名字重复的原因再决定是否需要按一下方式都保留

    > rownames_data <- make.names(data[,1],unique=T)  ##有效名称由字母 数字及dot和下划线组成,字符串中的空格或者-等“非法”字符都会被转为.。
    > length(rownames_data)  ##上面讲第一列导出到rownames_data,遇上重复名称,按顺序加上.1 .2 ...等后缀。
    [1] 27186

    > data <- data[,-1,drop = F]  ##[1,]访问第一行or 1:n;[,1]访问第一列or 1:n。这里是得到去掉第一列的剩余数据,drop=F表示更改是余下数据类型不变。

    > rownames(data) <- rownames_data  ##将没有重复的第一列数据保存到data中,作为行名。

    ****

    > write.table(rownames_data,file = "rownames_data.csv",quote = F)  ##不加quote参数,字符串及行号都会加上引号。行号默认有1 2 3。。列名x
    > write.table(rownames_data,file = "rownames_data_1.csv",quote = F,row.names = F)  ##
    > write.table(rownames_data,file = "rownames_data_1.csv",quote = F,row.names = F,col.names = F)  ##不要引号,不要行列名

    ****

    > data <- data[rowSums(data)>0,] #去掉和为零的行

    > nrow(data)
    [1] 27186


    > data <- data[apply(data, 1,var)!= 0,] ##!=表示不等于;去掉方差为0 的行,这些本身没有意义,也妨碍后续运算;1表示对阵列的行进行操作

    > ?rowSums  ##计算一个阵列的行的和,返回的是一个向量,长度为这个阵列的行数。
    > x <- cbind(x1 = 3, x2 = c(4:1, 2:5))
    > x
    x1 x2
    [1,] 3 4
    [2,] 3 3
    [3,] 3 2
    [4,] 3 1
    [5,] 3 2
    [6,] 3 3
    [7,] 3 4
    [8,] 3 5
    > rowSums(x)
    [1] 7 6 5 4 5 6 7 8
    > colSums(x)
    x1 x2
    24 24
    > x1 <- x[rowSums(x)>5,]
    > x1
    x1 x2
    [1,] 3 4
    [2,] 3 3
    [3,] 3 3
    [4,] 3 4
    [5,] 3 5

    ## apply(X, MARGIN, FUN, ...) 函数:Returns a vector or array or list of values obtained by applying a function to margins of an array or matrix.

    ##

    X
    an array, including a matrix.

    MARGIN
    a vector giving the subscripts which the function will be applied over. E.g., for a matrix 1 indicates rows, 2 indicates columns, c(1, 2) indicates rows and columns. Where X has named dimnames, it can be a character vector selecting dimension names.

    ****

    > mads <- apply(data,1,mad)  ##计算中值绝对偏差 (MAD, median absolute deviation)度量基因表达变化幅度

    > data <- data[rev(order(mads)),] ##将数据先按索引排序-升序order(),然后在反过来,然后在按照反过来的索引顺序形成data。

    > a <- c(1,4,0,9,57,32,6)
    > order(a)
    [1] 3 1 2 7 4 6 5  ##0在a中的索引值是3,故3在order()之后的第一位,order()返回的是索引值,以此类推
    
    > rev(order(a))
    [1] 5 6 4 7 2 1 3
    > b <- a[rev(order(a))]
    > b <- a[rev(order(a)),]
    Error in a[rev(order(a)), ] : 量度数目不对
    > b <- a[rev(order(a))]
    > b
    [1] 57 32 9 6 4 1 0
    > abc <- data.frame(a = c(1:3),b = c(4:6),c = c(7:9))
    > abc
      a b c
    1 1 4 7
    2 2 5 8
    3 3 6 9
    > mads <-  apply(abc,1,mad)
    > mads
    [1] 4.4478 4.4478 4.4478
    > cbind(abc,d = c(1,1,1))
      a b c d
    1 1 4 7 1
    2 2 5 8 1
    3 3 6 9 1
    > abc
      a b c
    1 1 4 7
    2 2 5 8
    3 3 6 9
    > abc <-  cbind(abc,d = c(1,1,1)) ##给加一列
    > abc
      a b c d
    1 1 4 7 1
    2 2 5 8 1
    3 3 6 9 1
    > mads <- apply(abc,1,mad)
    > abc
      a b c d
    1 1 4 7 1
    2 2 5 8 1
    3 3 6 9 1
    > mads
    [1] 2.2239 2.9652 3.7065
    > abc  <- abc[rev(order(mads)),] ##从结果看,按行的mad的降序重新排列了数据框
    > abc
      a b c d
    3 3 6 9 1
    2 2 5 8 1
    1 1 4 7 1

    ##data[x:y],访问数据框时,如果不加标点访问的是啥呢?

    > data_t[1:5]
    [1] 245667.7 427435.1 221687.5 371144.2 240187.2
    > data_t[,1:5]
                       FN1      DCN     CEMIP    CCDC80    IGFBP5
    untrt_N61311  245667.7 212953.1  40996.34 137229.15  77812.65
    untrt_N052611 427435.1 360796.2 137783.10 232772.17 288609.20
    untrt_N080611 221687.5 258977.3  53813.92  86258.13 210628.87
    untrt_N061011 371144.2 408573.1  91066.80 212237.32 168067.42
    trt_N61311    240187.2 210002.2  62301.12 136730.76  96021.74
    trt_N052611   450103.2 316009.1 223111.85 226070.89 217439.21
    trt_N080611   280226.2 225547.4 212724.84 124634.56 162677.38
    trt_N061011   376518.2 393843.7 157919.47 236237.81 168387.36

    > data_t <- t(data)  ##t()函数对数据框进行转置,行变列+列变行,对角线不变。

    > abc_t <- t(abc)
    > abc
      a b c d
    3 3 6 9 1
    2 2 5 8 1
    1 1 4 7 1
    > abc_t
      3 2 1
    a 3 2 1
    b 6 5 4
    c 9 8 7
    d 1 1 1

    ##文中说rows are samples and columns are variables,这句话不太懂。

    > variableL <- ncol(data_t)
    > variableL
    [1] 27186

    > if(sampleFile != ""){  ##这里的if()是何意???
    +     sample <- read.table(sampleFile,header = T,row.names = 1,sep ="	")  
    +     data_t_m <- merge(data_t,sample,by = 0)
    +     rownames(data_t_m) <- data_t_m$Row.names
    +     data_t <- data_t_m[,-1]
    + }
    > data_t[,1:5]
                       FN1      DCN     CEMIP    CCDC80    IGFBP5
    trt_N052611   450103.2 316009.1 223111.85 226070.89 217439.21
    trt_N061011   376518.2 393843.7 157919.47 236237.81 168387.36
    trt_N080611   280226.2 225547.4 212724.84 124634.56 162677.38
    trt_N61311    240187.2 210002.2  62301.12 136730.76  96021.74
    untrt_N052611 427435.1 360796.2 137783.10 232772.17 288609.20
    untrt_N061011 371144.2 408573.1  91066.80 212237.32 168067.42
    untrt_N080611 221687.5 258977.3  53813.92  86258.13 210628.87
    untrt_N61311  245667.7 212953.1  40996.34 137229.15  77812.65

    ## row.names = 1?

    ## merge()?

    ## rownames(data_t_m) <- data_t_m$Row.names?

    row.names = 1?

    > sample_2 <- read.table("new_1.txt",header = T,row.names = 2,sep ="	")
    > sample_2
                    Samp
    untrt1  untrt_N61311
    untrt2 untrt_N052611
    untrt3 untrt_N080611
    untrt4 untrt_N061011
    trt5      trt_N61311
    trt6     trt_N052611
    trt7     trt_N080611
    trt8     trt_N061011
    > sample_2 <- read.table("new_1.txt",header = T,row.names = NULL,sep ="	")
    > sample_2
               Samp conditions
    1  untrt_N61311     untrt1
    2 untrt_N052611     untrt2
    3 untrt_N080611     untrt3
    4 untrt_N061011     untrt4
    5    trt_N61311       trt5
    6   trt_N052611       trt6
    7   trt_N080611       trt7
    8   trt_N061011       trt8
    ##看来参数row.names是指定哪一列作为行名。
    ##所以上面的row.names = 1表示的是以第一列为行名,所以原数据框中第一行第一列的Samp就去掉了。
    
    

    merge()?

    > data_t[,1:5]
                       FN1      DCN     CEMIP    CCDC80    IGFBP5
    untrt_N61311  245667.7 212953.1  40996.34 137229.15  77812.65
    untrt_N052611 427435.1 360796.2 137783.10 232772.17 288609.20
    untrt_N080611 221687.5 258977.3  53813.92  86258.13 210628.87
    untrt_N061011 371144.2 408573.1  91066.80 212237.32 168067.42
    trt_N61311    240187.2 210002.2  62301.12 136730.76  96021.74
    trt_N052611   450103.2 316009.1 223111.85 226070.89 217439.21
    trt_N080611   280226.2 225547.4 212724.84 124634.56 162677.38
    trt_N061011   376518.2 393843.7 157919.47 236237.81 168387.36
    > sample <- read.table(sampleFile,header = T,row.names = 1,sep ="	")
    > sample
                  conditions
    untrt_N61311       untrt
    untrt_N052611      untrt
    untrt_N080611      untrt
    untrt_N061011      untrt
    trt_N61311           trt
    trt_N052611          trt
    trt_N080611          trt
    trt_N061011          trt ##这里可以看出data_t与sample的行名是一致的。隐隐约约觉得有个问题:这俩只是行名一致,又不是某列,这也能合并啊??
    > data_t_m <- merge(data_t,sample,by = 0)  ##merge()中参数all默认为FALSE,即合并时取交集--返回匹配的行;all = TRUE时,取按相同列取并集。
    > data_t_m[,1:5]  ##data_t_m的前5列    ###为何相同的列合并之后有了行名,还叫Row.names??
          Row.names      FN1      DCN     CEMIP    CCDC80
    1   trt_N052611 450103.2 316009.1 223111.85 226070.89
    2   trt_N061011 376518.2 393843.7 157919.47 236237.81
    3   trt_N080611 280226.2 225547.4 212724.84 124634.56
    4    trt_N61311 240187.2 210002.2  62301.12 136730.76
    5 untrt_N052611 427435.1 360796.2 137783.10 232772.17
    6 untrt_N061011 371144.2 408573.1  91066.80 212237.32
    7 untrt_N080611 221687.5 258977.3  53813.92  86258.13
    8  untrt_N61311 245667.7 212953.1  40996.34 137229.15
    > ncol(sample)
    [1] 1
    > ncol(data_t)
    [1] 27186
    > ncol(data_t_m)
    [1] 27188
    > data_t_m[,27182:27188] ##后7列
    FSTL4 COL23A1 BEST2 IBSP LTF TAC1 conditions
    1 1.515631 0.0000000 0.000000 0.0000000 1.722211 0.000000 trt 2 0.000000 1.7616370 1.181002 0.0000000 0.000000 0.000000 trt 3 0.000000 0.7821581 1.935371 0.7084139 2.398343 0.000000 trt 4 0.000000 3.7734772 2.274773 0.0000000 0.000000 3.956904 trt 5 5.018991 0.0000000 0.000000 4.1753179 0.000000 0.960596 untrt 6 0.000000 0.0000000 0.000000 0.0000000 0.000000 0.000000 untrt 7 2.512400 0.0000000 0.000000 0.0000000 2.629531 0.000000 untrt 8 0.000000 0.0000000 0.000000 0.0000000 0.000000 0.000000 untrt

    ####不明白的地方是data_t和sample合并之后为何行的顺序会发生变化,且变化之后的顺序又是什么顺序呢?

    rownames(data_t_m) <- data_t_m$Row.names?  

    > rownames(data_t_m) <- data_t_m$Row.names  ##把Row.names这一列(是第一列)内容设置为行名,但是这一列不变!
    > data_t <- data_t_m[,-1]  ##去掉第一列,

    > pca <- prcomp(data_t[,1:variableL],scale = T)

    PCA结果展示  ##有限展示,因为刚接触,这部分图的逻辑不太懂

    > fviz_eig(pca, addlabels = TRUE) ##addlabel参数为T,显示百分比。

    > fviz_pca_ind(pca, repel=T) ##repel=T 自动调整文本位置
    > fviz_pca_ind(pca)

    后面略过,没跟着学了

    参考https://mp.weixin.qq.com/s/4R14xJkQVPtaufaoXOcPIw

  • 相关阅读:
    [转帖]ExtJs与服务器的交互(一)
    Ext_两种处理服务器端返回值的方式
    命令行调用dubbo远程服务
    【转】Lombok Pojo默认初始值问题
    【转】Python之可变参数,*参数,**参数,以及传入*参数,进行解包
    Python获取并输出当前日期时间
    idea上 实现了Serializable接口,要自动生成serialVersionUID的方法
    mac下python2.x和python3.x的安装方法和升级方法/卸载
    工程代码不编译src的java目录下的xml文件问题及解决
    Spring Aop织入点语法
  • 原文地址:https://www.cnblogs.com/SWTwanzhu/p/12932858.html
Copyright © 2020-2023  润新知