生信在线作图-学习生物信息学有哪些比较好的网站或论坛？

ATAC-seq信息分析流程主要分为以下几个部分:数据质控、序列比对、峰检测、motif分析、峰注释、富集分析,下面将对各部分内容进行展开讲解。

下机数据经过过滤去除接头含量过高或低质量的reads,得到clean reads用于后续分析。常见的trim软件有Trimmomatic、Skewer、fastp等。fastp是一款比较新的软件,使用时可以用--adapter_sequence/--adapter_sequence_r2参数传入接头序列,也可以不填这两个参数,软件会自动识别接头并进行剪切。如:

fastp \

--in1 A1_1. \ # read1原始fq文件

--out1 A1_clean_1. \ # read1过滤后输出的fq文件

--in2 A1_2.  \ # read2原始fq文件

--out2 A1_clean_2. \ # read2过滤后输出的fq文件

--cut_tail  \ #从3’端向5’端滑窗,如果窗口内碱基的平均质量值小于设定阈值,则剪切

--cut_tail_window_size=1 \ #窗口大小

--cut_tail_mean_quality=30 \ #cut_tail参数对应的平均质量阈值

--average_qual=30 \ #如果一条read的碱基平均质量值小于该值即会被舍弃

--length_required=20  \ #经过剪切后的reads长度如果小于该值会被舍弃

fastp软件的详细使用方法可参考:

。fastp软件对于trim结果会生成网页版的报告,可参考官网示例

和

,也可以用FastQC软件对trim前后的数据质量进行评估,FastQC软件会对单端的数据给出结果,如果是PE测序需要分别运行两次来评估read1和read2的数据质量。

如:

fastqc A1_1.

fastqc A1_2.

FastQC会对reads从碱基质量、接头含量、N含量、高重复序列等多个方面对reads质量进行评估,生成详细的网页版报告,可参考官网示例:

经过trim得到的reads可以使用BWA、bowtie2等软件进行比对。首先需要确定参考基因组fa文件,对fa文件建立索引。不同的软件有各自建立索引的命令,BWA软件可以参考如下方式建立索引:

bwa index

建立好索引后即可开始比对,ATAC-seq推荐使用mem算法,输出文件经samtools排序输出bam:

bwa mem  A1_clean_1. A1_clean_2.

| samtools sort -O bam -T A1 >

值得注意的是,在实验过程中质体并不能完全去除,因此会有部分reads比对到质体序列上,需要去除比对到质体上的序列,去除质体序列可以通过samtools提取,具体方法如下:首先将不含质体的染色体名称写到一个chrlist文件中,一条染色体的名称写成一行,然后执行如下命令即可得到去除质体的bam

samtools view -b $chrlist > _MT_

用于后续分析的reads需要时唯一比对且去重复的,bwa比对结果可以通过MAPQ值来提取唯一比对reads,可以用picard、sambamba等软件去除dup,最终得到唯一比对且去重复的bam文件。

比对后得到的bam文件可以转化为bigWig(bw)格式,通过可视化软件进行展示。deeptools软件可以实现bw格式转化和可视化展示。首先需要在linux环境中安装deeptools软件,可以用以下命令实现bam向bw格式的转换:

bamCoverage -b -o

此外,可以使用deeptools软件展示reads在特定区域的分布,如:

computeMatrix reference-point   \ # reference-pioint表示计算一个参照点附近的reads分布,与之相对的是scale-regions,计算一个区域附近的reads分布

--referencePoint TSS   \#以输入的bed文件的起始位置作为参照点

-S   \ #可以是一个或多个bw文件

-R   \ #基因组位置文件

-b 3000   \ #计算边界为参考点上游3000bp

-a 3000   \ #计算边界为参考点下游3000bp,与-b合起来就是绘制参考点上下游3000bp以内的reads分布

-o   \ #输出作图数据名称

#图形绘制

plotHeatmap \

-m  new_ \ #上一步生成的作图数据

-out \ # 输出图片名称

绘图结果展示:

MACS2能够检测DNA片断的富集区域,是ATAC-seq数据call peak的主流软件。峰检出的原理如下:首先将所有的reads都向3'方向延伸插入片段长度,然后将基因组进行滑窗,计算该窗口的dynamic λ,λ的计算公式为:λlocal = λBG(λBG是指背景区域上的reads数目),然后利用泊松分布模型的公式计算该窗口的显著性P值,最后对每一个窗口的显著性P值进行FDR校正。默认校正后的P值(即qvalue)小于或者等于0.05的区域为peak区域。需要现在linux环境中安装macs2软件,然后执行以下命令:

macs2 callpeak \

-t \ #bam文件

-n A1 \ # 输出文件前缀名

--shift -100 \ #extsize的一半乘以-1

--extsize 200 \ #一般是核小体大小

--call-summits #检测峰顶信息

注:以上参数参考文献(Jie Wang,.2018.“ATAC-Seq analysis reveals a widespread decrease of chromatin accessibility in age-related macular degeneration.”Nature Communications)

ATAC分析得到的peak是染色质上的开放区域,这些染色质开放区域常常预示着转录因子的结合,因此对peak区域进行motif分析很有意义。常见的motif分析软件有homer和MEME。以homer软件为例,首先在linux环境中安装homer,然后用以下命令进行motif分析:

\

A1_ \ #用于进行motif分析的bed文件

 \ #参考基因组fa文件

A1  \ #输出文件前缀

-size  given \ #使用给定的bed区域位置进行分析,如果填-size -100,50则是用给定bed中间位置的上游100bp到下游50bp的区域进行分析

homer分析motif的原理及结果参见:

根据motif与已知转录因子的富集情况可以绘制气泡图,从而可以看到样本与已知转录因子的富集显著性。

差异peak代表着比较组合染色质开放性有差异的位点,ChIP-seq和ATAC-seq都可以用DiffBind进行差异分析。DiffBind通过可以通过bam文件和peak的bed文件计算出peak区域标准化的readcount,可以选择edgeR、DESeq2等模型进行差异分析。

在科研分析中我们往往需要将peak区域与基因联系起来,也就是通过对peak进行注释找到peak相关基因。常见的peak注释软件有ChIPseeker、homer、PeakAnnotator等。以ChIPseeker为例,需要在R中安装ChIPseeker包和GenomicFeatures包,然后就可以进行分析了。

library(ChIPseeker)

library(GenomicFeatures)

txdb<- makeTxDbFromGFF(‘’)#生成txdb对象,如果研究物种没有已知的TxDb,可以用GenomicFeatures中的函数生成

peakfile <-readPeakFile(‘A1_’)#导入需要注释的peak文件

peakAnno <- annotatePeak(peakfile,tssRegion=c(-2000, 2000), TxDb=txdb)

# 用peak文件和txdb进行peak注释,这里可以通过tssRegion定义TSS区域的区间

对于peak注释的结果,也可以进行可视化展示,如:

p <- plotAnnoPie(peakAnno)

通过注释得到的peak相关基因可以使用goseq、topGO等R包进行GO富集分析,用kobas进行kegg富集分析,也可以使用DAVID在线工具来完成富集分析。可以通过挑选感兴趣的GO term或pathway进一步筛选候选基因。

主要分为三个步骤:数据获取、数据处理、数据分析。数据获取是通过各种方式从实验中获得数据的过程。数据处理是对获得的数据进行预处理、标准化、格式化等步骤的过程。而数据分析则是对处理后的数据进行统计、图形展示、差异分析、关联分析等。

要想成功入门生信分析,首先要学会使用相应的工具来完成上述步骤。常用的工具包括R语言和Python语言。R语言是一门专门用于统计分布和图形展示的语言,而Python语言则通用性强,能够用于各个不同步骤中。学习R语言和Python语法时可以通过在线教程或者书本来学习。

在学习了R语言和Python语法之后就可以使用相应的工具来安装生信相关的包了。常用的包有bioconductor和anaconda。生信分析是一门非常庞大的学科,如果你想要入门,可以考虑从学习素材开始。书籍方面,《生信实战》是一本很好的入门书,而如果你想要更加深入学习生信分析,可以考虑阅读《生信数据分析》。在线课程方面,“生信进阶之路”是我个人非常喜欢的一个免费课程,它包含了很多实用的生信分析方法。另外还有“基因组学”、“微生物学”、“遗传学”这三门课程,都可以帮助你更好地理解生信分析。

生信菜鸟团上大学之后,我上网找资料时发现的第一个博客就是生信菜鸟团,里面包罗万象,涵盖很多方面(初次发现时,就感觉自己进入了新的天地)rabbit gao's blog  我超喜欢这个师兄的博客里面的笔记,很直观,尤其是python那部分。他是以代码的形式展示内容。沈梦圆博客梦圆师姐,和我一样喜欢用熊猫头像,她的博客也是刚刚建立不长时间。师姐的文笔很赞,看里面博文相信对你有帮助的。生信日志|鸣一道鸣一道师兄的博客我比较喜欢的是R做图那一块plob这个我比较少看,不过内容也不错,我后续再写上这个博客的描述。陈连福博客 听说连福老师有开培训班,实力自然也不差。糗世界←欢迎来到糗糗的世界糗世界主要包括:序列比对与NGS  R/bioconductorcircos教程,其中糗世界关于R和bioconductor以及NGS的归纳总结特别详尽生信客部落生信客部落是我自己的博客,刚建不久(2016.9.3建的),我目前在准备考研,打理的时间不多。但相信是一只潜力股,有提升的空间。也欢迎博友们交换 "友情链接".hope博客  hope 他(她)有一篇关于生物信息学在线工具的总结,我特别喜欢科研动力“endnote使用宝典”,专注写endnote相关的内容。(注:endnote 是文献管理的软件,插入引用文献的神器)biochen生物伯臣生物里也蛮多归纳整理的Bob's Blog   bob这位兄弟的博客我接触不多,我后续补上描述.论坛(包括生信论坛和其他一些相关的网站):生信技能树生信技能树前面那个师兄有详细描述过。我也亲眼见证了它从无到有的过程,看着生信技能树感觉特别亲切,感觉就像自家的孩子一样。我自己由于准备考研和书写毕业论文的事情,在生信技能树建设的参与度不高。总之,好喜欢......生物信息学天空内容超全的一个生信论坛丁香园(生信板块)丁香园,就不解释了,一个国内最成功的论坛之一。医学生基本都知道的一个论坛。小木虫小木虫,里面蛮多资源的,也是国内最成功的论坛之一生物统计家园描述待输入...基因堂描述待输入biostars这是一个生信问答网站生信刷题网站ROSALIND | About 这个是一个生物信息的刷题网站,超多实战题(纯英文,既提高英语水平,又训练了自己的实战能力,何乐而不为)。实战走起......生物信息学在线工具网站生信客部落生物信息工具整合(包含在线工具与离线工具)C更新中这里面包含了一些生物信息学可视化工具,包含在线的工具和一些离线的可视化工具,由于目前个人水平有限,所以还有待继续完善