纳米孔Nanopore-16S数据分析学习笔记

最新推荐文章于 2023-10-20 09:21:17 发布

原创最新推荐文章于 2023-10-20 09:21:17 发布 · 2.9k 阅读

6 ·

CC 4.0 BY-SA版权

文章标签：

#纳米孔 16s

生物信息同时被 2 个专栏收录

59 篇文章

订阅专栏

测序

4 篇文章

订阅专栏

本文详细介绍了纳米孔测序数据从下载、预处理到分析的全过程，包括使用INC-seq流程获得共识序列，以及共识序列的方向校正和OTU表的获取。过程中解决了数据下载和格式转换的问题，并探讨了不同比对方法的效果。

摘要生成于 C知道，由 DeepSeek-R1 满血版支持，前往体验 >

2020.5找到了开放的protocol：https://www.protocols.io/view/nanoampli-seq-bioinformatics-workflow-u25eyg6

已经转换成fasta的数据供练手：http://gigadb.org/dataset/100530

1.下载原始数据

本次学习分析的文章是这篇：https://academic.oup.com/gigascience/article/7/12/giy140/5202451
这篇文章的原始数据有点问题，使用sra和ena数据库直接下载都基本上会失败，sra只能下到一个10M左右的数据，转换格式成fastq后只能获得4.6M的数据。最后使用aspera connect下载可以成功。命令如下，我是黑果，其他系统格式应该类似，软件安装和使用参见我前面的学习记录。https://jiawen.zd200572.com/916.html

ascp -QT -l 3000m -P33001 -L- -k1 -i asperaweb_id_dsa.openssh era-fasp@fasp.sra.ebi.ac.uk:vol1/run/ERR224/ERR2241542/RAW_R9.5_1_organism.tar.gz BioData/Nanopore/

下载的数据没有通过完整性测试，没办法进行下一步分析。后面重新下载了原始数据，见下面：

2.解压转换成fastq

ftp://parrot.genomics.cn/gigadb/pub/10.5524/100001_101000/100530/nanopore-master.zip
#解压
unzip nanopore-master.zip
#hd5转换成fastq

1.把下载的fastq格式转为fasta

其实这是一个很简单的过程，即使自己用个脚本或者使用命令行也能解决，介于想要重复作者结果，就按作者的原步骤进行。

#新建一个python2的conda环境
conda create -n python2 python=2
#安装依赖
conda install  -y  biopython=1.65 blast seqtk numpy vsearch mafft-ginsi
# seqtk转换fastq到fasta
seqtk seq  ERR2241540.fastq -a > ERR2241540.fasta

2.获得共识序列

这里走了点弯路，其实本文的参考文献里说明了是使用INC-seq这个流程进行前处理的，找到这个流程的github仓库，就可以使用了。

git clone https://github.com/CSB5/INC-Seq.git
############我是分界线，以下是我走的弯路#######
#安装三代纠错软件，卡在这好久，docker下载会卡住失败，源码编译会下载一个不存在的boost_1_58_0-headersonly.tbz2文件，最后用conda解决（友情提示，只有linux-64版本）
conda install pbdagcon -c https://mirrors.tuna.tsinghua.edu.cn/anaconda/cloud/bioconda/
#########这些分界线包着的是不需要的步骤，都包含在了INC-seq的脚本里####
#开始获得共识序列#
./inc-seq.py -i  ERR2241540.fasta -o consensus.fa  --copy_num_thre 3  --iterative  -a poa
#INC-seq的流程包含了三个比对方法（blastn,  graphmap, poa），默认是采用blastn的，依照文章作者的参数，--copy_num_thre 3  --iterative ，只有poa方法能够获得共识序列，blastn提示：
    #Warning: [blastn] Examining 5 or more matches is recommended
    #Warning: [blastn] Examining 5 or more matches is recommended
    #Warning: [blastn] Examining 5 or more matches is recommended
    #Warning: [blastn] Examining 5 or more matches is recommended
    #PBDAGCON failed (trimmed more than 100 bases)!
    #Consensus construction failed!
#graphmap提示Candidate read found!
    #PBDAGCON failed (trimmed more than 100 bases)!
    #Consensus construction failed!
#完成后的文件列表在这里，当然，只有poa有结果。可以看出序列利用率是比较低的，好多序列由于长度不够，或者发现的片段不一致而过滤掉了。
tree -h 
.
├── [ 183]  bpipe.config
├── [   0]  consensus-blastn.fa
├── [   0]  consensus-gra.fa
├── [150K]  consensus-poa.fasta
├── [4.0K]  data
│   └── [ 17K]  inc_seq_test_read.fa
├── [2.3M]  ERR2241540.fasta
├── [4.6M]  ERR2241540.fastq
├── [7.4K]  inc-seq.py
├── [1.4K]  LICENSE
├── [1.2K]  pipeline.bpipe
├── [3.0K]  README.md
└── [4.0K]  utils

3.共识序列方向校正

#最近发现需要安装几个依赖，在脚本文件里写得比较清楚了，这里补上
wget  ftp://emboss.open-bio.org/pub/EMBOSS/emboss-latest.tar.gz
tar zxvf emboss-latest.tar.gz
cd EMBOSS-6.6.0
./configure
make -j 5
make install
wget https://drive5.com/python/python_scripts.tar.gz
tar zxvf python_scripts.tar.gz
sudo mv uc2otutab.py /bin
sudo chmod +x /bin/uc2otutab.py
chopSEQ.py -i inc_seq_test_read.fa -f <forward_primer> -r <reverse_primer> > <filtered_fasta_file>

4.获得otu表

awk -v k="Sample1" '/^>/{gsub(">","",$0);$0=">barcodelabel="k";"$0}1' chop.fasta > Sample1_barcode.fa

python nanoCLUST.py -i Sample1_barcode.fa -o ./nano -s 0,450,451,900,901,1300

比较遗憾，尽管没有报错，却没有得出结果，otu表为空，或许是原始数据有问题，或许前面参数出了问题，做了个探索。

#2020.5更新，添加依赖，应该可以正常运行出结果了。