導(dǎo)讀
需求:1. 讀取一個(gè)目錄下的所有數(shù)據(jù)表;2. 統(tǒng)計(jì)特定元素的非"NA"的頻數(shù),并匯總到一個(gè)表中匾浪。方法:使用R語(yǔ)言read.table讀取數(shù)據(jù)表昼伴,list保存數(shù)據(jù)表噩咪,sum統(tǒng)計(jì)元素頻數(shù)弄慰,write.table保存結(jié)果序芦、for循環(huán)批量操作每一個(gè)數(shù)據(jù)表等等。
一蔬咬、查看數(shù)據(jù)表
Linux中查看文件信息
ls -alh
# 查看所有文件名,如下:
總用量 3.9M
drwxrwxr-x 2 cheng WST 4.0K 10月 18 15:46 .
drwxrwxr-x 8 cheng WST 4.0K 10月 16 16:12 ..
-rw-rw-r-- 1 cheng WST 197K 10月 16 11:00 bin.25.tsv
-rw-rw-r-- 1 cheng WST 122K 10月 16 11:00 bin.33.tsv
-rw-rw-r-- 1 cheng WST 176K 10月 16 11:00 bin.36.tsv
-rw-rw-r-- 1 cheng WST 153K 10月 16 11:00 bin.37.tsv
-rw-rw-r-- 1 cheng WST 260K 10月 16 11:00 bin.38.tsv
-rw-rw-r-- 1 cheng WST 253K 10月 16 11:00 bin.39.tsv
-rw-rw-r-- 1 cheng WST 217K 10月 16 11:00 bin.43.tsv
-rw-rw-r-- 1 cheng WST 81K 10月 16 11:00 bin.46.tsv
-rw-rw-r-- 1 cheng WST 228K 10月 16 11:00 bin.50.tsv
-rw-rw-r-- 1 cheng WST 114K 10月 16 11:00 bin.56.tsv
-rw-rw-r-- 1 cheng WST 237K 10月 16 11:00 bin.58.tsv
-rw-rw-r-- 1 cheng WST 152K 10月 16 11:00 bin.61.tsv
-rw-rw-r-- 1 cheng WST 238K 10月 16 11:00 bin.63.tsv
-rw-rw-r-- 1 cheng WST 144K 10月 16 11:00 bin.65.tsv
-rw-rw-r-- 1 cheng WST 102K 10月 16 11:00 bin.66.tsv
-rw-rw-r-- 1 cheng WST 371K 10月 16 11:00 bin.67.tsv
-rw-rw-r-- 1 cheng WST 87K 10月 16 11:00 bin.71.tsv
-rw-rw-r-- 1 cheng WST 198K 10月 16 11:00 bin.81.tsv
-rw-rw-r-- 1 cheng WST 114K 10月 16 11:00 bin.89.tsv
-rw-rw-r-- 1 cheng WST 171K 10月 16 11:00 bin.91.tsv
-rw-rw-r-- 1 cheng WST 279K 10月 16 11:00 bin.94.tsv
head bin.25.tsv
# 查看文件格式沐寺,如下:
locus_tag ftype length_bp gene EC_number COG product
LBILEGMC_00001 CDS 324 hypothetical protein
LBILEGMC_00002 CDS 2589 tmoS_1 2.7.13.3 Sensor histidine kinase TmoS
LBILEGMC_00003 CDS 852 hypothetical protein
LBILEGMC_00004 CDS 1164 hypothetical protein
LBILEGMC_00005 CDS 1356 hypothetical protein
LBILEGMC_00006 CDS 975 2.5.1.10 COG0142 (2E,6E)-farnesyl diphosphate synthase
LBILEGMC_00007 CDS 708 hypothetical protein
LBILEGMC_00008 CDS 780 dtd3_1 3.1.1.96 D-aminoacyl-tRNA deacylase
LBILEGMC_00009 CDS 699 guaA_1 6.3.5.2 GMP synthase [glutamine-hydrolyzing]
二林艘、讀取文件名
list.files:讀取工作目錄下所有文件名
strsplit:切割字符串
as.character:將strsplit結(jié)果【list】轉(zhuǎn)成字符【character】
for循環(huán):批處理
files=list.files(pattern="bin.*.tsv")
# 讀取所有文件的全名,含后綴
name=vector()
for(i in 1:length(files))
{
name[i]=as.character(strsplit(files[i], split=".tsv"))
# 批量提取文件名
}
name
# 查看文件名混坞,如下:
[1] "bin.25" "bin.33" "bin.36" "bin.37" "bin.38" "bin.39" "bin.43" "bin.46"
[9] "bin.50" "bin.56" "bin.58" "bin.61" "bin.63" "bin.65" "bin.66" "bin.67"
[17] "bin.71" "bin.81" "bin.89" "bin.91" "bin.94"
length(name)
# 一共21個(gè)文件
[1] 21
三狐援、讀取數(shù)據(jù)表/框
for循環(huán):批處理
list:用列表保存所有數(shù)據(jù)表信息
read.table:讀表函數(shù)
sep='\t':制表符分隔【參數(shù)】
na.string="":空為NA【參數(shù)】
stringsAsFactors=F:不把字符串的列辨認(rèn)成因子【參數(shù)】
header=T:第一行為表頭【參數(shù)】
quote="":無(wú)引號(hào)【參數(shù)】
comment.char="":無(wú)注釋【參數(shù)】
ml=list()
# 創(chuàng)建一個(gè)新列表啥酱,存放每個(gè)數(shù)據(jù)框的數(shù)據(jù)。
for(i in 1:length(files))
{
ml[[i]]=read.table(files[i], sep='\t', na.string="", stringsAsFactors=F, header=T, quote="", comment.char="")
# 讀取所有數(shù)據(jù)框
}
length(ml)
# 一共21個(gè)
[1] 21
head(ml[[1]])
# 挑一個(gè)檢查厨诸,沒(méi)啥毛病镶殷,如下:
locus_tag ftype length_bp gene EC_number COG
1 LBILEGMC_00001 CDS 324 <NA> <NA> <NA>
2 LBILEGMC_00002 CDS 2589 tmoS_1 2.7.13.3 <NA>
3 LBILEGMC_00003 CDS 852 <NA> <NA> <NA>
4 LBILEGMC_00004 CDS 1164 <NA> <NA> <NA>
5 LBILEGMC_00005 CDS 1356 <NA> <NA> <NA>
6 LBILEGMC_00006 CDS 975 <NA> 2.5.1.10 COG0142
product
1 hypothetical protein
2 Sensor histidine kinase TmoS
3 hypothetical protein
4 hypothetical protein
5 hypothetical protein
6 (2E,6E)-farnesyl diphosphate synthase
四、統(tǒng)計(jì)“gene”非"NA"的數(shù)量
head(ml[[1]]["gene"])
# 查看需要統(tǒng)計(jì)的信息微酬,計(jì)算非<NA>的基因數(shù)量:
gene
1 <NA>
2 tmoS_1
3 <NA>
4 <NA>
5 <NA>
6 <NA>
gene_num=vector()
# 創(chuàng)建一個(gè)新向量绘趋,存放列表中每個(gè)數(shù)據(jù)框的gene數(shù)量。
for(i in 1:length(files))
{
gene_num[i]=lengths(ml[[i]]["gene"])-sum(is.na(ml[[i]]["gene"]))
# 基因總數(shù)/BIN:
# is.na判斷是否為NA
# gene數(shù)/數(shù)據(jù)框=gene列長(zhǎng)度/數(shù)據(jù)框-NA數(shù)/數(shù)據(jù)框
}
length(gene_num)
# 每個(gè)數(shù)據(jù)框一個(gè)gene數(shù)
[1] 21
gene_num
# 查看每個(gè)數(shù)據(jù)框的“gene”數(shù)量颗管,如下:
[1] 1559 1303 1755 1319 2054 2636 2372 806 1972 1023 1754 1445 2276 1287 896
[16] 3349 839 1593 1016 1750 2397
五陷遮、結(jié)果合并、保存
data.frame:可合并多個(gè)向量成表
write.table:保存數(shù)據(jù)表/框
sep="\t":制表符分隔【參數(shù)】
row.names=F:去行標(biāo)【參數(shù)】
quote=F:去引號(hào)【參數(shù)】
name
# 查看:
[1] "bin.25" "bin.33" "bin.36" "bin.37" "bin.38" "bin.39" "bin.43" "bin.46"
[9] "bin.50" "bin.56" "bin.58" "bin.61" "bin.63" "bin.65" "bin.66" "bin.67"
[17] "bin.71" "bin.81" "bin.89" "bin.91" "bin.94"
gene_num
# 查看:
[1] 1559 1303 1755 1319 2054 2636 2372 806 1972 1023 1754 1445 2276 1287 896
[16] 3349 839 1593 1016 1750 2397
result=data.frame(name, gene_num)
# 合并
result
# 查看:
name gene_num
1 bin.25 1559
2 bin.33 1303
3 bin.36 1755
4 bin.37 1319
5 bin.38 2054
6 bin.39 2636
7 bin.43 2372
8 bin.46 806
9 bin.50 1972
10 bin.56 1023
11 bin.58 1754
12 bin.61 1445
13 bin.63 2276
14 bin.65 1287
15 bin.66 896
16 bin.67 3349
17 bin.71 839
18 bin.81 1593
19 bin.89 1016
20 bin.91 1750
21 bin.94 2397
write.table(result, file="gene.tsv", sep="\t", row.names=F, quote=F)
# 保存結(jié)果到當(dāng)前工作目錄
Linux中查看gene.tsv文件
cat gene.tsv
# 結(jié)果沒(méi)毛病垦江,如下:
name gene_num
bin.25 1559
bin.33 1303
bin.36 1755
bin.37 1319
bin.38 2054
bin.39 2636
bin.43 2372
bin.46 806
bin.50 1972
bin.56 1023
bin.58 1754
bin.61 1445
bin.63 2276
bin.65 1287
bin.66 896
bin.67 3349
bin.71 839
bin.81 1593
bin.89 1016
bin.91 1750
bin.94 2397