首先做完了周末班工作, 包括软件安装以及R包安装:
打开 Rstudio告诉我它的工作目录。
getwd()
新建6个向量,基于不同的原子类型
。(重点是字符串,数值,逻辑值)
a<-c('March','April','May')
b<-c(1,7,3,5)
c<-c(F,T,T)
d<-c(1:21)
e<-paste0(rep("sample",3),1:3)
告诉我在你打开的rstudio里面 getwd()
代码运行后返回的是什么?
新建一些数据结构,比如矩阵,数组,数据框,列表等。重点是数据框,矩阵)
1)矩阵
m<-matrix(1:21,ncol =7)
m
rownames(m)<-paste0(rep("gene",3),1:3)
colnames(m)<-paste0(rep("sample",7),1:7)
m
矩阵:向量+维度属性,这个维度是一个整数向量,只包含两个元素:行数和列数,分别用nrow和ncol表示。
1)用matrix()函数创建矩阵
2)用矩阵维度的属性dim()查看矩阵的行列数。
3)用rbind()函数按行拼接矩阵; 用cbind()函数按列拼接矩阵
2)数组
array
数组与矩阵类似,但是数组的维度可以大于2。
1、用array()函数创建2维数组
ar <- array(1:24,dim=c(4,6)
这就产生了一个4行6列的二维数组,也是4行6列的矩阵。
2、用array()函数创建3维数组
ar <- array(1:18,c(2,3,3))
3)数据框
df<-data.frame(gene=paste0("gene",1:5),a1 = rnorm(n=5,mean=1,sd=0.5),a2 = rnorm(n=5,mean=2),a3 = rnorm(n=5,mean=3),a4 = rnorm(n=5,mean=4),a5 = rnorm(n=5,mean=5))
df
在你新建的数据框进行切片操作,比如首先取第1,3行, 然后取第4,6列
df[c(1,3),] ##行
df[,c(4,6)] ##列
4)列表
mylist <- list(stud.id=1234,stud.name = "Tom", stud.marks = c(12,3,14,25,19))
mylis
使用data函数来加载R内置数据集 rivers
, 描述它。
data(rivers)
rivers
下载 https://www.ncbi.nlm.nih.gov/sra?term=SRP133642 里面的 RunInfo Table
文件读入到R里面,了解这个数据框,多少列,每一列都是什么属性的元素。(参考B站生信小技巧获取runinfo table) 这是一个单细胞转录组项目的数据,共768个细胞,如果你找不到RunInfo Table
文件,可以点击下载,然后读入你的R里面也可以。
sra<-read.table(file="SraRunTable.txt",header=T,sep=',')
View(sra)
dim(sra)
ncol(sra)
class(sra[,1])
class(sra[1,])
mode(sra[,2])
for (i in 1:ncol(sra)) {
x=class(sra[,i])
print(x)
}
下载 https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE111229 里面的样本信息sample.csv
读入到R里面,了解这个数据框,多少列,每一列都是什么属性的元素。(参考 https://mp.weixin.qq.com/s/fbHMNXOdwiQX5BAlci8brA 获取样本信息sample.csv)如果你实在是找不到样本信息文件sample.csv,也可以点击下载。
sample<-read.csv(file="sample.csv")
##也可以用 sample<- read.table("sample.csv",sep = ",")
View(sample)
colnames(sample)
dim(sample)
ncol(sample)
for (i in 1:ncol(sample)){
y=class(sample[,i])
print (y)
}
把前面两个步骤的两个表(RunInfo Table
文件,样本信息sample.csv)关联起来,使用merge函数。
dim(sra)
dim(sample)
merge <- merge(sra,sample,by.x = "Sample.Name",by.y = "Accession")
# merge() 函数 此时sra的Sample.Name列和sample的Accession列内容相同,合并这一列
dim(merge)
View(merge)
基于下午的统计可视化
对前面读取的 RunInfo Table
文件在R里面探索其MBases列,包括 箱线图(boxplot)和五分位数(fivenum),还有频数图(hist),以及密度图(density) 。
boxplot(sra$MBases, main = "boxplot")
plot(fivenum(sra$MBases), main = "fivenum")
hist(sra$MBases, main = "hist")
plot(density(sra$MBases,na.rm=T), main = "density")
把前面读取的样本信息
表格的样本名字根据下划线分割
看第3列元素的统计情况。第三列代表该样本所在的plate
title=sample$Title
class(title)
title=as.character(sample$Title)
class(title)
#修改title的格式为character
strsplit(title[1],split = "_")
# strsplit( )函数用于字符串分割,其中split 是分割参数。所得结果以默认以list形式展示。
plate <- unlist(lapply(title,function(x){strsplit(x,'_')[[1]][3]}))
## lapply 代表list 将指定操作应用于列表中的所有元素
# 这里写了一个自编函数,对每个title值进行上一行的操作,获得768个值
【没搞懂这个自编函数,关键是lapply出来的是list格式,用循环函数搞不定这个】
table(plate)
# table() 函数各个值出现的频率
根据plate把关联到的 RunInfo Table
信息的MBases列分组检验是否有统计学显著的差异。
t.test(merge$MBases~plate)
##666,要记住这个函数和格式。
分组绘制箱线图(boxplot),频数图(hist),以及密度图(density) 。
boxplot(merge$MBases~plate,main='boxplot')
e<-merge[,c("MBases","Title")]
e$plate = plate
View (e)
hist(e$MBases,main='boxplot',freq = F, breaks = "sturges")
plot(density(e$MBases))
##并不能实现分组绘制啊。。。。
使用ggplot2把上面的图进行重新绘制。
一个ggplot2教程
ggplot(data = mpg, aes(x = cty, y = hwy)) + ### 设置x, y以及数据
geom_point(aes(color = cyl)) + ### 添加 点, 以及 设置 点 的颜色
geom_smooth(method = "lm") + ### 添加拟合曲线/曲线以及置信区间
coord_cartesian() + ### 笛卡尔坐标系, 即平面直角坐标系
scale_color_gradient() + ### 图例的颜色
theme_base() + ### 主题为base风格
scale_shape(solid = FALSE) + ### 点的类型为中空。
ggtitle("New Plot Title") + ### 标题
xlab("New X label") + ### 横轴
ylab("New Y label") ### 纵轴
library(ggplot2)
#ggplot2将所绘制图形的各部分独立出来,如xy坐标,基本组件的颜色、图标大小、填充类型、堆叠方式、坐标轴、地图投影、数据分组、图形分面等等信息划归为一些基本类型,每一部分分别用一段代码表示。每段代码内部有相应的参数控制,各代码段再通过“+”运算符连接。这里的“+”并非四则运算的加法,而是表示图形各要素/组件之间的叠加。“+”运算符所连接的代码片段,有任何一部分做出更改,整个图形就做出相应调整。因此,本质上,ggplot2代码是用R语言实现的一种绘图语言
##会改代码就得了
e$num <- 1:768
colnames(e)
ggplot(e,aes(x=plate,y=MBases)) + geom_boxplot()
# 箱图
ggplot(e,aes(x=MBases)) + geom_histogram() + facet_grid(plate ~ .)
ggplot(e,aes(x=MBases,fill=plate))+geom_histogram(colour="grey")
# 频数图
ggplot(e,aes(x=MBases)) + geom_density() + facet_grid(plate ~ .)
ggplot(e,aes(x=MBases,fill=plate))+geom_density()
# 密度图
使用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"))
##太难了,会改就得了