一个引号引发的血案（文本文件读取的小事故）

生信技能树

发布于 2024-11-21 09:06:34

3200

安排学徒探索了一下表达量芯片的不同探针平台信息，然后学徒给我反馈了一个在他看来有意思的bug，就是在读取一个txt文件的时候会出现读不完整的情况：

k = read.table('./GPL570-55999.txt',
               header = T,
               sep = '\t',
               fill = T,
               #quote = '',
               check.names = F,
               skip = 16) 

Warning message:
In scan(file = file, what = what, sep = sep, quote = quote, dec = dec,  :
  EOF within quoted string

有时候Warning 信息也得注意

这个时候函数给出来的是Warning message，并不是大家常见的error，所以有些人可能会忽略它。

其实解决方式很简单，需要仔细调整参数，比如在上面的read.table 函数里面添加了 quote = '' 的参数，或者使用高级的R包：

k2 = rio::import('./GPL570-55999.txt',skip = 16)
k1 = data.table::fread('./GPL570-55999.txt',skip = 16)

上面的失败读取的问题在于最后一个探针 "AFFX-r2-Bs-dap-M_at" ：

> dim(k);dim(k1);dim(k2)
[1] 28059    16
[1] 54675    16
[1] 54675    16

> head(k[,1]);tail(k[,1])
[1] "1007_s_at" "1053_at"   "117_at"    "121_at"    "1255_g_at"
[6] "1294_at"  
[1] "AFFX-M27830_M_at"    "AFFX-PheX-3_at"      "AFFX-PheX-5_at"     
[4] "AFFX-PheX-M_at"      "AFFX-r2-Bs-dap-3_at" "AFFX-r2-Bs-dap-M_at"

所以需要去看看这个探针对应的注释信息，这个时候不用R代码了，就一个简简单单shell命令即可：

grep -n  'AFFX-r2-Bs-dap-M_at' GPL570-55999.txt 
54668:AFFX-r2-Bs-dap-M_at       --Control       Homo sapiens    Oct 6, 2014 Control sequence Affymetrix Proprietary Database B. subtilis /GEN=dapB, jojF /DB_XREF=gb:L38424.1 /NOTE=SIF corresponding to nucleotides 2055-2578 of gb:L38424.1 /DEF=Bacillus subtilis dihydropicolinate reductase (jojE) gene, complete cds; poly(A) polymerase (jojI) gene, complete cds; biotin acetyl-CoA-carboxylase ligase (birA) gene, complete cds; jojC, jojD, jojF, jojG, jojH genes, complete cds's. AFFX-r2-Bs-dap-M
-ne

其实我本来是准备让搜索如何输出指定的第 28059 行附近的信息，还使用人工智能大模型帮我搜索了方法（在Linux的shell里面如何针对txt文件输出指定的某一行）：

在Linux的shell中，可以使用几种不同的方法来输出文本文件中的指定行。以下是一些常用的命令：

使用sed命令：
- sed 'Nq' filename.txt：输出文件filename.txt的第N行。
使用awk命令：
- awk 'NR==N {print}' filename.txt：输出文件filename.txt的第N行，其中NR是awk内置的记录数变量。
使用head和tail命令（当你知道行号在文件的前几行时）：
- head -n N filename.txt | tail -n 1：输出文件filename.txt的第N行，适用于行号较小的情况。
使用grep命令（如果行包含特定的文本）：
- grep 'pattern' filename.txt：输出包含特定文本pattern的所有行。使用-n选项可以显示行号。
使用less或more命令（交互式查看）：
- 可以先使用less或more命令查看文件内容，然后使用上下箭头键跳转到指定的行。
使用nl命令（带行号查看）：
- nl filename.txt：显示文件filename.txt的所有行，并带上行号。然后你可以手动查找到指定的行。
使用cat与-n选项：
- cat -n filename.txt：类似于nl，cat -n命令也会输出文件的每一行，并在行首显示行号。
使用tac命令（从文件末尾开始计数）：
- tac filename.txt | grep '^[N]'：这将从文件末尾开始计数，显示倒数第N行。注意，^在正则表达式中表示行的开始，[N]是十进制数。

选择哪种方法取决于你的具体需求和文件的大小。对于大型文件，建议使用awk或sed，因为它们在处理大文件时更高效。对于较小的文件，或者当你只需要查看文件的一部分时，head、tail或grep可能更简单快捷。

首先，这个大模型给出来的第一个方法就是错误的，然后呢我自己提出来的输出指定的第 28059 行附近的信息这个思路是错误的，因为read.table 函数读取txt后的数据框的大小并不是txt文件本身出问题的地方。因为我们选择在上面的read.table 函数里面添加了 quote = '' 的参数，其实就明白了是因为txt文件里面有单引号，这个才是引发问题的原因，文本文件里面的有 1333 个地方都有这个单引号：

grep "'" GPL570-55999.txt |wc
    1333

有些时候连Warning都没有但也是错误的

比如另外一个GPL的txt文件读取：

k2 = rio::import('./GPL13158-5065.txt',skip = 16)
k1 = data.table::fread('./GPL13158-5065.txt',skip = 16)
k = read.table('./GPL13158-5065.txt',
               header = T,
               sep = '\t',
               fill = T,
               check.names = F,
               skip = 16)

完全是顺利的运行了：

> dim(k);dim(k1);dim(k2)
[1] 26374    16
[1] 54715    16
[1] 54715    16
> head(k[,1]);tail(k[,1])
[1] "1007_PM_s_at" "1053_PM_at"   "117_PM_at"    "121_PM_at"   
[5] "1255_PM_g_at" "1294_PM_at"  
[1] "AFFX-ThrX-3_at"  "AFFX-ThrX-5_at"  "AFFX-ThrX-M_at" 
[4] "AFFX-TrpnX-3_at" "AFFX-TrpnX-5_at" "AFFX-TrpnX-M_at"
>

高级函数也会报错

前面提到了用高级的R包：

k2 = rio::import('./GPL570-55999.txt',skip = 16)
k1 = data.table::fread('./GPL570-55999.txt',skip = 16)

但是它们也不是万能的哦，之前就遇到了一个单细胞转录组表达量矩阵文件，是txt或者csv格式的，使用fread就只能读取一半的基因或者细胞。

本文参与腾讯云自媒体同步曝光计划，分享自微信公众号。

原始发表：2024-08-01，如有侵权请联系 cloudcommunity@tencent.com 删除

搜索

本文分享自生信技能树微信公众号，前往查看

如有侵权，请联系 cloudcommunity@tencent.com 删除。

本文参与腾讯云自媒体同步曝光计划，欢迎热爱写作的你一起参与！

登录后参与评论

0 条评论

热度