在小丫畫(huà)圖交付的一個(gè)代碼項(xiàng)目中榨乎，需要先從XENA下載一個(gè)表達(dá)量數(shù)據(jù):https://toil.xenahubs.net/download/tcga_RSEM_gene_tpm.gz
樣本大概是10,5,35個(gè), 考慮到人類(lèi)的基因大概有2w多個(gè)，那么這就是一個(gè)10000 X 20000的大樣本數(shù)據(jù)辆琅，鑒于這還是一個(gè)TPM蕴茴，數(shù)據(jù)類(lèi)型是浮點(diǎn)型迅诬，文件解壓縮之后就是4.61G, 如果全部加載到R語(yǔ)言中箭养，大部分的電腦估計(jì)都受不了
library(pryr)
test <- data.table::fread("./tcga_RSEM_gene_tpm.gz")
object_size(test)
# 5.11 GB
考慮到并非所有數(shù)據(jù)都是我們所需要的奥喻，是否可以只讀取部分的數(shù)據(jù)呢偶宫？原作者的解決方案是通過(guò)R調(diào)用命令行的方式，提取部分?jǐn)?shù)據(jù)环鲤，然后讓R語(yǔ)言進(jìn)行加載纯趋。
可是大部分人的操作系統(tǒng)都是Windows，所有運(yùn)行的時(shí)候就會(huì)報(bào)錯(cuò)冷离，能不能就用戶(hù)R語(yǔ)言解決這個(gè)問(wèn)題呢吵冒？當(dāng)然可以，只要你認(rèn)真讀過(guò)read.table的那么多參數(shù)西剥，你就會(huì)知道他的那么多參數(shù)并不是裝飾用的痹栖。
讀取前幾行
讓我們先學(xué)習(xí)一個(gè)簡(jiǎn)單的參數(shù)nrows, 他的作用就是讀取前N行，知道它之后瞭空，那就不需要去調(diào)用head了
headtcga <- read.table("./tcga_RSEM_gene_tpm",
sep = "\t",
stringsAsFactors = FALSE,
nrow = 1)
效果就是讀取第一行揪阿，構(gòu)建一個(gè)數(shù)據(jù)框，然后將其轉(zhuǎn)成向量咆畏。但既然目標(biāo)是向量南捂，其實(shí)還有另一種實(shí)現(xiàn)方案，readLines讀取的就是一個(gè)字符串旧找，然后將其分隔成向量即可溺健。
headtcga <- readLines("tcga_RSEM_gene_tpm", n =1)
headtcga <- strsplit(headtcga, split="\t")[[1]]
讀取指定列
讀取指定列會(huì)稍微困難一些，因?yàn)?code>colClasses不太好理解钦讳。R語(yǔ)言在用read.table讀取數(shù)據(jù)的時(shí)候其實(shí)做了很多事情矿瘦，有一件事情就是負(fù)責(zé)確認(rèn)每一列的數(shù)據(jù)類(lèi)型枕面，R語(yǔ)言需要根據(jù)不同數(shù)據(jù)類(lèi)型進(jìn)行內(nèi)存分配。
如果你想實(shí)現(xiàn)讀取指定列缚去，那么你就得自己去設(shè)置每一列的數(shù)據(jù)類(lèi)型潮秘。如果哪些列不需要，就將其它的數(shù)據(jù)類(lèi)型定義為NULL易结，R語(yǔ)言就會(huì)忽略它枕荞。
讀取代碼如下：
cat(paste0("Begin at ", Sys.time(),"\n"))
first_5_rows <- read.table("./tcga_RSEM_gene_tpm", nrows = 5,
stringsAsFactors = FALSE,
header = FALSE,
skip = 1,
check.names = FALSE)
classes <- sapply(first_5_rows, class)
# targetnum 你需要讀取的列
classes[-targetnum] <- rep("NULL", length(classes) - length(targetnum)) #將非目標(biāo)列定義為NULL
classes[1] <- "character" # 加上第一列
# 讀取文件(跳過(guò)第一行)
targetCancerTPM <- read.table("tcga_RSEM_gene_tpm",
sep= "\t",
skip = 1,
colClasses = classes)
colnames(targetCancerTPM) <- tcgasample[targetnum]
targetCancerTPM[1:3, 1:3]
cat(paste0("End at ", Sys.time(),"\n"))
如果僅讀取我們需要的列的話(huà)，最終只消耗了500M的內(nèi)存搞动，相對(duì)于之前的5G內(nèi)存躏精，減少了將近10倍。
讀取指定行和指定列
這就是需要對(duì)文件進(jìn)行逐行讀取解析了鹦肿，我用readLines造了一個(gè)輪子矗烛，函數(shù)名為read_part,目前能用的參數(shù)為
- file: 輸入的文件路徑，支持.gz文件
- rows: 讀取指定行, 比如說(shuō)1:100, 就是前100行箩溃。當(dāng)為-1時(shí)則是讀取所有行
- rows: 讀取指定列, 比如說(shuō)
c(1,3,4,5,6), 就是1,3,4,5,6列瞭吃。當(dāng)為-1時(shí)則是讀取所有列 - comment.char = "#", 會(huì)把"#"開(kāi)頭的行忽略掉，這個(gè)參數(shù)我還需要考慮下是否保留涣旨。
# 函數(shù)目標(biāo):
# 讀取文件中的指定行和指定列
# 不包括注釋行
read_part <- function(file, rows = 1, columns = -1, sep = "\t",
stringsAsFactors = FALSE,
header = FALSE,
check.names = FALSE,
comment.char = "#", ...){
dfl <- list()
if (grepl("gz$", file)){
con <- gzfile(file, open = "rb")
} else{
con <- file(file, open = "r")
}
i <- 0
j <- 1
repeat{
rec <- readLines(con, 1)
if (length(rec) == 0) break
i <- i + 1
# 當(dāng)rows = -1時(shí), 會(huì)讀取所有行
# 超過(guò)目標(biāo)行時(shí)停止讀取
if (i > max(rows) & rows != -1) break
# 不考慮注釋行
if (grepl(comment.char, rec )) next
if ( ! i %in% rows & rows != -1) next
items <- strsplit(rec, split = sep, fixed = TRUE)[[1]]
if ( columns == -1){
select_cols <- items
} else{
select_cols <- items[columns]
}
#print(select_cols)
dfl[[j]] <- select_cols
j <- j + 1
}
close(con)
df <- do.call(rbind, dfl)
return(df)
}
