bwa bowtie2 salmon subread hisat2建索引和比对

bwa

# 建索引
bwa index -a bwtsw ~/hg38.fa

# 100bp长的双端测序reads的比对用bwa的mem算法
bwa mem -t 5 -R "@RG\tID:$sample\tSM:$sample\tLB:WGS\tPL:Illumina"  /public/biosoft/GATK/resources/bundle/hg38/bwa_index/gatk_hg38 7E5241.L1_1.fastq 7E5241.L1_2.fastq   | samtools sort -@ 5 -o 7E5241.bam -
cat config |while read id
do
    arr=($id) 
    fq1=${arr[1]} 
    fq2=${arr[2]} 
    sample=${arr[0]} 
    bwa mem -t 5 -R "@RG\tID:$sample\tSM:$sample\tLB:WGS\tPL:Illumina" $INDEX $fq1 $fq2  | samtools sort -@ 5 -o $sample.bam -   
done

bowtie

# 官网下载索引
# 自己构建
bowtie2-build   /home/qmcui/database/reference/hg38/hg38.fa   hg38 1>hg38.bowtie2.log   2>hg38.bowtie2.error   &
# Total time for backward call to driver() for mirror index: 01:18:29

# 增加线程数
bowtie2-build   /home/qmcui/database/reference/hg38/hg38.fa   hg38   --threads   5 
# Total time for backward call to driver() for mirror index: 00:21:55

# 单端比对
$bin_bowtie2  -p 5  -x  $bowtie2_index -U  $id | samtools sort  -O bam  -@ 5 -o - > ${sample}.bam 
# 双端序列比对
bowtie2  -p 6  -x  /public/reference/index/bowtie/hg38 -1 /home/qmcui/project/7.CHIP/2.val.fq/K_K4Me_clean_1_val_1.fq.gz -2 /home/qmcui/project/7.CHIP/2.val.fq/K_K4Me_clean_2_val_2.fq.gz | samtools sort  -O bam -@ 6 -o - >K_K4Me.sort.bam

salmon

直接从fastq到达特异性表达矩阵

# wget ftp://ftp.ensembl.org/pub/release-91/fasta/homo_sapiens/cdna/Homo_sapiens.GRCh38.cdna.all.fa.gz
# 注意,其他软件都是利用参考基因组的fa文件建立索引,这个软件用的cdna的数据建立索引
# 建索引
 ./bin/salmon index -t transcripts.fa -i transcripts_index
# 4min结束
# 双端比对
salmon quant -i transcripts_index -l <LIBTYPE> -1 reads1.fq -2 reads2.fq -o transcripts_quant
# 单端比对
salmon quant -i transcripts_index -l <LIBTYPE> -r reads.fq -o transcripts_quant
# 参考:https://salmon.readthedocs.io/en/latest/salmon.html#using-salmon

subread

# 建索引
subread-build index ../../genome/hg38/hg38.fa -o /home/qmcui/database/index/subread/hg38

subjunc -T 2  -i /home/qmcui/database/index/subread/hg38 -r ${id}_1_val_1.fq -R ${id}_2_val_2.fq -o ${id}.subjunc.sam 

备注:需要对应建索引的软件版本,版本不同会报错!

hisat2

# hisat官网
http://ccb.jhu.edu/software/hisat2/index.shtml
# 2018.11.30
# 直接官网下载,不用自己构建
wget -c ftp://ftp.ccb.jhu.edu/pub/infphilo/hisat2/data/hg38.tar.gz

# 自己建立
hisat2-build  –p  4  hg38.fa  hg38

# 比对
# hisat2 [options]* -x <hisat2-idx> {-1 <m1> -2 <m2> | -U <r> | –sra-acc <SRA accession number>} [-S <hit>]
# https://ccb.jhu.edu/software/hisat2/manual.shtml
hisat2 -p 16 -x ./genome_tran -1 SRR534293_1.fastq -2 SRR534293_2.fastq –S SRR534293.sam

HISAT2建立索引时,就应该把转录组信息加进去。
HISAT2提供两个Python脚本将GTF文件转换成hisat2-build能使用的文件:

extract_exons.py Homo_sapiens.GRCh38.83.chr.gtf > genome.exon
extract_splice_sites.py Homo_sapiens.GRCh38.83.chr.gtf > genome.ss

此外,HISAT2还支持将SNP信息加入到索引中,这样比对的时候就可以考虑SNP的情况。这仍然需要将SNP文件转换成hisat2-build能使用的文件:

extract_snps.py snp142Common.txt > genome.snp

最后,将基因组、转录组、SNP建立索引:

hisat2-build -p 4 genome.fa --snp genome.snp --ss genome.ss --exon genome.exon genome_snp_tran
官网下载展示

压缩包有make_grch38_tran.sh文件,详细记录了创建索引的过程。

$ cat make_hg38.sh 
#!/bin/sh

#
# Downloads sequence for the HG38 version of H. spiens (human) from
# UCSC.
#
# The base files, named ??.fa.gz
#
# By default, this script builds and index for just the base files,
# since alignments to those sequences are the most useful.  To change
# which categories are built by this script, edit the CHRS_TO_INDEX
# variable below.
#

UCSC_HG38_BASE=http://hgdownload.cse.ucsc.edu/goldenPath/hg38/bigZips
F=hg38.chromFa.tar.gz

get() {
    file=$1
    if ! wget --version >/dev/null 2>/dev/null ; then
        if ! curl --version >/dev/null 2>/dev/null ; then
            echo "Please install wget or curl somewhere in your PATH"
            exit 1
        fi
        curl -o `basename $1` $1
        return $?
    else
        wget $1
        return $?
    fi
}

HISAT2_BUILD_EXE=./hisat2-build
if [ ! -x "$HISAT2_BUILD_EXE" ] ; then
    if ! which hisat2-build ; then
        echo "Could not find hisat2-build in current directory or in PATH"
        exit 1
    else
        HISAT2_BUILD_EXE=`which hisat2-build`
    fi
fi

rm -f genome.fa
get ${UCSC_HG38_BASE}/$F || (echo "Error getting $F" && exit 1)
tar xvzfO $F > genome.fa || (echo "Error unzipping $F" && exit 1)
rm $F

CMD="${HISAT2_BUILD_EXE} genome.fa genome"
echo Running $CMD
if $CMD ; then
    echo "genome index built; you may remove fasta files"
else
    echo "Index building failed; see error message"
fi
hisat2下载解压

salmon

软件版本older,也会报错!跟subread一样。

# 查看可用的命令
# Salmon v0.9.1
salmon -h
# 查看帮助文档之Salmon's quasi-mapping-based mode
salmon --no-version-check quant --help-reads
# 查看帮助文档之Salmon's alignment-based mode
salmon --no-version-check quant --help-alignment
wget ftp://ftp.ensembl.org/pub/release-91/fasta/homo_sapiens/cdna/Homo_sapiens.GRCh38.cdna.all.fa.gz
wget ftp://ftp.ensembl.org/pub/release-91/fasta/homo_sapiens/cdna/Homo_sapiens.GRCh38.cdna.abinitio.fa.gz
wget ftp://ftp.ensembl.org/pub/release-91/fasta/mus_musculus/cdna/Mus_musculus.GRCm38.cdna.all.fa.gz
wget ftp://ftp.ensembl.org/pub/release-91/fasta/mus_musculus/cdna/Mus_musculus.GRCm38.cdna.abinitio.fa.gz
salmon index -t /home/jianmingzeng/reference/transcriptome/Homo_sapiens.GRCh38.cdna.all.fa.gz  -i hg38_index
salmon index -t /home/jianmingzeng/reference/transcriptome/Mus_musculus.GRCm38.cdna.all.fa.gz   -i mm10_index
最后编辑于
©著作权归作者所有,转载或内容合作请联系作者
  • 序言:七十年代末,一起剥皮案震惊了整个滨河市,随后出现的几起案子,更是在滨河造成了极大的恐慌,老刑警刘岩,带你破解...
    沈念sama阅读 204,590评论 6 478
  • 序言:滨河连续发生了三起死亡事件,死亡现场离奇诡异,居然都是意外死亡,警方通过查阅死者的电脑和手机,发现死者居然都...
    沈念sama阅读 86,808评论 2 381
  • 文/潘晓璐 我一进店门,熙熙楼的掌柜王于贵愁眉苦脸地迎上来,“玉大人,你说我怎么就摊上这事。” “怎么了?”我有些...
    开封第一讲书人阅读 151,151评论 0 337
  • 文/不坏的土叔 我叫张陵,是天一观的道长。 经常有香客问我,道长,这世上最难降的妖魔是什么? 我笑而不...
    开封第一讲书人阅读 54,779评论 1 277
  • 正文 为了忘掉前任,我火速办了婚礼,结果婚礼上,老公的妹妹穿的比我还像新娘。我一直安慰自己,他们只是感情好,可当我...
    茶点故事阅读 63,773评论 5 367
  • 文/花漫 我一把揭开白布。 她就那样静静地躺着,像睡着了一般。 火红的嫁衣衬着肌肤如雪。 梳的纹丝不乱的头发上,一...
    开封第一讲书人阅读 48,656评论 1 281
  • 那天,我揣着相机与录音,去河边找鬼。 笑死,一个胖子当着我的面吹牛,可吹牛的内容都是我干的。 我是一名探鬼主播,决...
    沈念sama阅读 38,022评论 3 398
  • 文/苍兰香墨 我猛地睁开眼,长吁一口气:“原来是场噩梦啊……” “哼!你这毒妇竟也来了?” 一声冷哼从身侧响起,我...
    开封第一讲书人阅读 36,678评论 0 258
  • 序言:老挝万荣一对情侣失踪,失踪者是张志新(化名)和其女友刘颖,没想到半个月后,有当地人在树林里发现了一具尸体,经...
    沈念sama阅读 41,038评论 1 299
  • 正文 独居荒郊野岭守林人离奇死亡,尸身上长有42处带血的脓包…… 初始之章·张勋 以下内容为张勋视角 年9月15日...
    茶点故事阅读 35,659评论 2 321
  • 正文 我和宋清朗相恋三年,在试婚纱的时候发现自己被绿了。 大学时的朋友给我发了我未婚夫和他白月光在一起吃饭的照片。...
    茶点故事阅读 37,756评论 1 330
  • 序言:一个原本活蹦乱跳的男人离奇死亡,死状恐怖,灵堂内的尸体忽然破棺而出,到底是诈尸还是另有隐情,我是刑警宁泽,带...
    沈念sama阅读 33,411评论 4 321
  • 正文 年R本政府宣布,位于F岛的核电站,受9级特大地震影响,放射性物质发生泄漏。R本人自食恶果不足惜,却给世界环境...
    茶点故事阅读 39,005评论 3 307
  • 文/蒙蒙 一、第九天 我趴在偏房一处隐蔽的房顶上张望。 院中可真热闹,春花似锦、人声如沸。这庄子的主人今日做“春日...
    开封第一讲书人阅读 29,973评论 0 19
  • 文/苍兰香墨 我抬头看了看天上的太阳。三九已至,却和暖如春,着一层夹袄步出监牢的瞬间,已是汗流浃背。 一阵脚步声响...
    开封第一讲书人阅读 31,203评论 1 260
  • 我被黑心中介骗来泰国打工, 没想到刚下飞机就差点儿被人妖公主榨干…… 1. 我叫王不留,地道东北人。 一个月前我还...
    沈念sama阅读 45,053评论 2 350
  • 正文 我出身青楼,却偏偏与公主长得像,于是被迫代替她去往敌国和亲。 传闻我的和亲对象是个残疾皇子,可洞房花烛夜当晚...
    茶点故事阅读 42,495评论 2 343

推荐阅读更多精彩内容