用于记录R基本知识，勤学勤练习
grep
# a为一个data frame,取含有RNA-Seq的行
index <- grep("RNA-Seq", a$Assay_Type)
b <- a[index,]
读取文件，选择不读哪一行
##不读取带有！的行
data <- read.table("file",comment.char="!", sep="\t")
去掉数据框的某一列,添加行名, 列名
##删除第一列, b是data.frame
b<-b[,-1]
##添加行名，列名
rowname(b) <- b[,1]
colname(b) <- b[1,]
sort
##a为一data frame
##最小值
sort(a$num)[1] ; 等价min(a$num)
##最大值
sort(a$num, decreasing=T)[1] ; 等价max(a$num)
##five number, 最小值，下四分位，中位数，上四分位，最大值
fivenum(a$num)
##<5000的行
b <- a[a$num < 5000,]
strsplit 切割
#以2进行分割
b <- "adfdfd2jdfkd"
strsplit(b, "2")
adfdf, jdfkd
t检验
t.test(1:10, y = c(7:20))
t = -5.4349, df = 21.982, p-value = 1.855e-05
#p value <0.05, 则显著
paste, rep
##a, 连续粘贴4次
paste("a",1:5)
#结果：
a 1, a 2, a 3, a 4
##a, 连续粘贴4次,以“-”连接
paste("a",1:5, sep='-')
结果：
a-1, a-2, a-3, a-4
##复制a 4次
rep("a", 4)
结果：
a,a,a,a
apply 循环
##b 为一个data frame, 每一行求平均数
apply(b, 1, function(x){
mean(x)
})
其中 1 代表行，2 代表列
上面可以简写为：
apply(b, 1, mean)
# 获取超过5个元素>0 的行
b[apply(b ,1 function(x) sum(x>1) >5), ]
输出图片
# 比如输出PDF图片
pdf("test.pdf", width=12, height=10)
ggplot(data, aes(x,y)) +geom_line()
dev.off
quantile()
## 在R语言中取百分位比用quantile()函数
> data <- c(1,2,3,4,5,6,7,8,9,10)
> quantile(data,0.5)
50%
5.5
> quantile(data,c(0.25,0.75))
25% 75%
3.25 7.75
## 可以画阈值线，比如BSA画图时;有index 数据 data
p + geom_hline(yintercept =quantile(data, 0.95) ,color="blue",linetype=2,size-0.5)
which()
# which 可以筛选数据库中特有的行，比如有数据框data，筛选出delt大于0的行
zheng <- data[which(data$delt>0),]
melt()，dcast()
## melt（） 和 dcast（）为reshape2中的命令
## melt()： 可以将长的数据框变为短的
> head(new_data)
pos xx_mean xy_mean yy_mean
19 36001 1.4052965 0.8884643 1.0787320
20 38001 0.8538059 1.0241777 1.1718369
21 40001 1.7150511 1.5263786 1.4663535
23 44001 0.0000000 1.8001673 0.9257189
24 46001 0.2837331 1.5285861 1.5567648
25 48001 0.6670272 2.0511266 0.0000000
> head(melt(new_data,id.vars = "pos"))
pos variable value
1 36001 xx_mean 1.4052965
2 38001 xx_mean 0.8538059
3 40001 xx_mean 1.7150511
4 44001 xx_mean 0.0000000
5 46001 xx_mean 0.2837331
6 48001 xx_mean 0.6670272
以“pos”为不变的ID列，其它变量合并在variable列中，其对应的值为value列
## dcast()： 可以将短的数据框变为长的
> head(test)
pos variable value
1 36001 xx_mean 1.4052965
2 38001 xx_mean 0.8538059
3 40001 xx_mean 1.7150511
4 44001 xx_mean 0.0000000
5 46001 xx_mean 0.2837331
6 48001 xx_mean 0.6670272
> head(dcast(test,pos~variable))
pos xx_mean xy_mean yy_mean
1 36001 1.4052965 0.8884643 1.0787320
2 38001 0.8538059 1.0241777 1.1718369
3 40001 1.7150511 1.5263786 1.4663535
4 44001 0.0000000 1.8001673 0.9257189
5 46001 0.2837331 1.5285861 1.5567648
6 48001 0.6670272 2.0511266 0.0000000
#～左边的表示不变的列，右边则展开
dplyr包处理数据
均以如下例子进行操作
> head(data)
Chrom Pos Pvalue
1 7 5813753 4.191874e-10
2 7 5831698 2.779123e-09
3 7 5818992 3.800454e-09
4 7 5850455 1.097069e-05
5 3 14961172 1.521800e-05
6 7 5808739 1.833734e-05
- filter()
筛选Chrom >6，Pos >14961172
data_chr <- data%>% filter(Chrom>6,Pos>14961172)
> head(data_chr)
Chrom Pos Pvalue
1 7 5813753 4.191874e-10
2 7 5831698 2.779123e-09
3 7 5818992 3.800454e-09
4 7 5850455 1.097069e-05
- select()
该函数可以去除某一列,去除Chrom那一列
> data_min <- data%>% select(-Chrom)
> head(data_min)
Pos Pvalue
1 5813753 4.191874e-10
2 5831698 2.779123e-09
3 5818992 3.800454e-09
4 5850455 1.097069e-05
5 14961172 1.521800e-05
6 5808739 1.833734e-05
- filter_all(), filter_if(), filter_at()
以上函数均对数值型数值操作，若不是，则需要用selec()去除
去除所有属性小于6的行
data_p <- data%>% filter_all(all_vars(.<6))
筛选任意一属性>3的行
data_p <- data%>% filter_all(any_vars(.>3))
筛选以seq开头的属性任一大于3的行
data_r <- data%>% filter_at(vars(starts_with("Seq")),any_vars(.>3))
筛选以d开头的属性任一可被2整除的行
data_R <- data%>% filter_at(vars(starts_with("Seq")),any_vars((. %%2)==0))
- group_by() summaries()
group_by()指在对一数据框按照一定变量进行分组，一般和summaries() 一起用
# 按cyl进行分组
by_cyl <- mtcars %>% group_by(cyl)
#
by_cyl %>% summarise(
disp = mean(disp),
hp = mean(hp)
)
# A tibble: 3 x 3
cyl disp hp
<dbl> <dbl> <dbl>
1 4 105. 82.6
2 6 183. 122.
3 8 353. 209
cyl组中有4，6，8三个变量，所以根据该变量进行求均值
也可以同时对两个变量进行分组
#现对vs进行分组，然后对am进行分组，vs中有0，1两个值，am中也有0，1两个值
by_vs_am <- mtcars %>% group_by(vs, am)
# 分别求数量
by_vs <- by_vs_am %>% summarise(n = n())
> by_vs
# A tibble: 4 x 3
# Groups: vs [2]
vs am n
<dbl> <dbl> <int>
1 0 0 12
2 0 1 6
3 1 0 7
4 1 1 7
# 可以看到vs 为0，am为0 共有12个，vs为0，am为1的有6个
mutata()
用于添加新的一列
> head(mtcars)
mpg cyl disp hp drat wt qsec vs am gear carb
Mazda RX4 21.0 6 160 110 3.90 2.620 16.46 0 1 4 4
Mazda RX4 Wag 21.0 6 160 110 3.90 2.875 17.02 0 1 4 4
Datsun 710 22.8 4 108 93 3.85 2.320 18.61 1 1 4 1
> head(mtcars %>% mutate(total =1))
mpg cyl disp hp drat wt qsec vs am gear carb total
1 21.0 6 160 110 3.90 2.620 16.46 0 1 4 4 1
2 21.0 6 160 110 3.90 2.875 17.02 0 1 4 4 1
3 22.8 4 108 93 3.85 2.320 18.61 1 1 4 1 1
cbind(),rbind()
合并两个数据框，类似于cat,paste
test <- cbind(data1,data2)
即可进行合并
rainbow()
得到不同的颜色
data <- read.table("IV_vs_non_IV",sep = "\t",header = T)
cols <- rainbow(5*1.2)
cols
[1] "#FF0000FF" "#FFFF00FF" "#00FF00FF" "#00FFFFFF" "#0000FFFF" "#FF00FFFF"
boxplot(ks ~ gene, data = data, col = cols)
更换两列内容
比如有一数据框
head(data)
Gene Go level
1 gene1 GO:1 BP
2 gene2 Go:2 MF
将第二列和第一列互换
head(data[c(2,1,3)])
Go Gene level
1 GO:1 gene1 BP
2 Go:2 gene2 MF
reorder() 排序
现有两列数据
v1 v2
1 23
2 32
3 2
4 9
## 以v2从大到小排序后，得到对应v1
reorder(v1,v2)
## 以v2从小到大排序后，得到对应v1
reorder(v1,-v2)
BiocManager() 可非常方便安装R包
首先安装BiocManager
install.packages('BiocManager')
BiocManager::install('')
R包pheatmap
- border_color
热图的单位元素的描边颜色，NA表示不描边
默认：“gray60” - cluster_rows (按原来顺序排列)
用布尔值来判定行元素是不是应该被clustered or clust object。
默认：TRUE
同理，cluster_cols决定是否应该把列clust
for 循环
library(ggplot2)
## 具体根据染色体数量进行调整；下面图一个3X4 共12条染色体
par(mfrow =c(3,4))
data <- read.table('Tri_ang_all_ltr.4buyi1Mb200kb_4r',sep = '\t',header = T)
对应染色体名称
vec <- c('LG01','LG02','LG03','LG04')
data1 <- data[which(data$Chr==val),]
barplot(P ~ BP, data = data1, col='#3366cc', border='#3366cc',ylab = 'Number',xlab = paste0(val,' (Mb)'))
}
查看R包路径，以及修改路径
1.1 更改R包安装路径
- 在R里面运行命令(假定修改为如下路径)
.libPaths('/data/pub//R/x86_64-redhat-linux-gnu-library/3.6')
- 在安装包是定义
install.packages('LDheapment', lib='/data/pub//R/x86_64-redhat-linux-gnu-library/3.6')
1.2 查看R包的路径
> .libPaths()
[1] "/data/pub//R/x86_64-redhat-linux-gnu-library/3.6"
[2] "/data/pub//miniconda2/lib/R/library"
可以看到有2个r包的路径
1.3 修改R包路径
- 临时修改
.libPaths('/data/pub//R/x86_64-redhat-linux-gnu-library/3.6')
临时讲上述路径添加为r包的路径
- 永久性修改
新建一个 .Rprofile的文本文件
touch .Rprofile
并将 要定义的路径存入该文本文档中。路径如下
.libPaths(c("path1","path2","path3"))
package ‘XXX’ is not available for this version of R
- R包开源网站
- 找到所需要安装的包，比如LDheatmap
- 安装
install.packages("https://cran.r-project.org/src/contrib/Archive/LDheatmap/LDheatmap_1.0-6.tar.gz",repos = NULL, type="source")
鉴定密度图的峰值
library(reshape2)
## 读取数据
data <- read.table('test.txt')
## 鉴定峰值函数
densFindPeak <- function(x){
td <- density(x)
maxDens <- which.max(td$y)
list(x=td$x[maxDens],y=td$y[maxDens])
}
#计算峰值
densFindPeak(data$v1)
不同的group分别统计fivenumer
## 有data数据
type value
A 34
A 343
A 45
B 43
B 34
B 65
C 7
C 9
C 10
## 统计
library(tapply)
tapply(df$data, df$groupBy, summary)
tapply(data$value, data$type, summary)
基本操作
packageVersion("rvcheck") # 查看rvcheck版本
## 删除R包
remove.packages("rvcheck")
# 下载包
download.file("https://cran.r-project.org/src/contrib/Archive/rvcheck/rvcheck_0.1.8.tar.gz","rvcheck_0.1.8.tar.gz")
## 本地安装包
install.packages("rvcheck_0.1.8.tar.gz",repos = NULL)
## 在线安装包
install.packages("https://cran.r-project.org/src/contrib/Archive/rvcheck/rvcheck_0.1.8.tar.gz",repos = NULL,type = "source")
ifelse()，将指定列数值大于一定数值进行替换
# 创建一个示例数据框
df <- data.frame(A = c(1,2,3,4,5), B = c(1,2,4,5,6))
# 使用ifelse函数，将第二列（B列）中大于3的值替换为3
df$B <- ifelse(df$B > 3, 3, df$B)
