打開 Rstudio 告訴我它的工作目錄
getwd()
新建6個(gè)向量,基于不同的原子類型免绿。(重點(diǎn)是字符串唧席,數(shù)值,邏輯值)
a <- c('apple','pear','hello')
b <- c(3,5,7,9)
c <- c(F,T)
告訴我在你打開的rstudio里面 getwd() 代碼運(yùn)行后返回的是什么嘲驾?
getwd()
新建一些數(shù)據(jù)結(jié)構(gòu)淌哟,比如矩陣,數(shù)組辽故,數(shù)據(jù)框徒仓,列表等重點(diǎn)是數(shù)據(jù)框,矩陣)
mat <- matrix(1:20,4,5)
array1 <- array(1:24,c(2,3,4))
patientid <- c(1,2,3,4)
age <- c(23,44,56,78)
etilogy <- c('CAD','MI','HF','CAD')
sur <- c(T,F,T,F)
EF <- c(25,45,22,40)
HF <- data.frame(patientid,age,etilogy,sur,EF)
myhead <-'My <- list'
g <- c(23,4,5,6)
h <- matrix(1:10,2)
k <- c("allen",'pear','watermelon')
mylist <- list(myhead=myhead,age_sg=g,h,k)
在你新建的數(shù)據(jù)框進(jìn)行切片操作誊垢,比如首先取第1掉弛,3行, 然后取第4喂走,6列
HF[1:3,]
HF[,4:5]
data(rivers)
rivers
下載 https://www.ncbi.nlm.nih.gov/sra?term=SRP133642 里面的 RunInfo Table 文件讀入到R里面,了解這個(gè)數(shù)據(jù)框芋肠,多少列乎芳,每一列都是什么屬性的元素。(參考B站生信小技巧獲取runinfo table) 這是一個(gè)單細(xì)胞轉(zhuǎn)錄組項(xiàng)目的數(shù)據(jù)帖池,共768個(gè)細(xì)胞奈惑,如果你找不到RunInfo Table 文件,可以點(diǎn)擊下載睡汹,然后讀入你的R里面也可以肴甸。
SRP <- read.table("SRP133642.txt",header = T,sep='\t')
dim(SRP)
class(SRP$BioSample)
mode(SRP$BioSample)
ncol(SRP)
1:ncol(SRP)
class(SRP[,1])
mode(SRP[,2])
for (i in 1:ncol(SRP)) {
x=class(SRP[,i])
print(x)
}
load("sample.Rdata")
for (i in 1:ncol(pdata)) {
x=class(pdata[,i])
print(x)
}
把前面兩個(gè)步驟的兩個(gè)表(RunInfo Table 文件焕檬,樣本信息sample.csv)關(guān)聯(lián)起來(lái), 使用merge函數(shù)澳泵。
m1 <- merge(pdata,a,by.y = 'Sample_Name',by.x = 'geo_accession')
基于下午的統(tǒng)計(jì)可視化实愚。對(duì)前面讀取的 RunInfo Table(SRP) 文件在R里面探索其MBases列,包括 箱線圖(boxplot)和五分位數(shù)(fivenum)兔辅,還有頻數(shù)圖(hist)腊敲,以及密度圖(density) 。
SRP$MBases
boxplot(SRP$MBases,main='boxplot')
plot(fivenum(SRP$MBases),main='fivenum')
hist(SRP$MBases,main = 'hist')
plot(density(SRP$MBases),main='density')
把前面讀取的樣本信息表格的樣本名字根據(jù)下劃線分割看第3列元素的統(tǒng)計(jì)情況维苔。第三列代表該樣本所在的plate
load(file="sample.Rdata")
Title=pdata$title
Title=as.character(pdata$title)
class(Title)
strsplit(Title[1],'_')[[1]][3]
length(Title)
plate <- unlist(lapply(Title,function(x){
x
strsplit(x,'_')[[1]][3]
}))
plate
根據(jù)plate把關(guān)聯(lián)到的 RunInfo Table 信息的MBases列分組檢驗(yàn)是否有統(tǒng)計(jì)學(xué)顯著的差異碰辅。
a <- read.table('SRP133642.txt',header = T,sep='\t')
save(a,file='RunInfo.Rdata')
m1 <- merge(pdata,a,by.y = 'Sample_Name',by.x = 'geo_accession')
t.test(m1$MBases~plate)
分組繪制箱線圖(boxplot),頻數(shù)圖(hist)介时,以及密度圖(density)
boxplot(m1$MBases~plate,main='boxplot')
e <- m1[,c("title",'MBases')]
e$plate <- plate
hist(e$MBases,main='boxplot',freq = F, breaks = "sturges")
plot(density(e$MBases))
使用ggplot2把上面的圖進(jìn)行重新繪制没宾。
library(ggplot2)
e$num <- 1:768
colnames(e)
ggplot(e,aes(x=plate,y=MBases)) + geom_boxplot()
# 頻數(shù)圖
ggplot(e,aes(x=MBases)) + geom_histogram(fill="lightblue",colour="grey") + facet_grid(plate ~ .) # 頻數(shù)圖
ggplot(e,aes(x=MBases,fill=plate))+geom_histogram()
ggplot(e,aes(x=MBases,fill=plate))+geom_histogram()
# 密度圖
ggplot(e,aes(y=MBases,x=num)) + geom_point() + stat_density2d(aes(alpha=..density..),geom = "raster",contour = F)+ facet_grid(plate ~ .)
ggplot(e,aes(x=MBases,fill=plate))+geom_density()
使用ggpubr把上面的圖進(jìn)行重新繪制。
# install.packages("ggpubr")
library(ggpubr)
ggboxplot(e, x="plate", y="MBases", color = "plate", palette = "aaas",add = "jitter") + stat_compare_means(method = "t.test")
gghistogram(e, x="MBases", fill = "plate",palette = c("#f4424e", "#41a6f4"))
ggdensity(e, x="MBases", fill = "plate", color = "plate", add = "mean",palette = c("#f4424e", "#41a6f4"))
隨機(jī)取384個(gè)MBases信息沸柔,跟前面的兩個(gè)plate的信息組合成新的數(shù)據(jù)框循衰,第一列是分組,第二列是MBases褐澎,總共是384*3行數(shù)據(jù)会钝。
sample(e$MBases,384)
nrow(e)
data <- e[sample(nrow(e),384),]
data <- data[,c(3,2,4,1)]
str(data)