18 thoughts on “ChIP-Seq 数据分析

  1. 欧老师,我现在得到chip-seq数据的一个表,就是Active region(就是不同样本之间的interval重合)的数据列表。实验一共是6个样本,其中3个为一类。
    列表的格式大概是
    第一列:Active region ID
    第二列: 染色体ID
    第三列: Active Region 开始的核苷酸位置
    第四列: Active Region 结束的核苷酸位置
    五到十列: 6个样本在Active region 内的峰值密度Peak fragment density
    请问现在如果我要寻找在两组样本里差异最大的Active region,是不是利用edgeR、DESeq2等方法处理峰值密度就可以得到结果?
    谢谢!

    1. 这里你需要注意的是Peak fragment density是什么样的数据,如果是raw reads counts的话,可以使用edgeR/DESeq2,但如果不是的话,你直接log2转换后使用limma、anova即可。

      1. 峰顶,Value at Summit (peak calling应该用的是MACS),所以可以直接用edgeR/DESeq2上?

        1. 理论上edgeR/DESeq2用的是峰宽所在区域的所有reads的总数,否则p-value会不准确。中间还小涉及到normalization的问题。我个人觉着还是不要使用为好。

          1. 欧老师,
            问题一
            请问我现在是只有两个数据 一个是均值
            Average of fragment densities of all the 32-nt length bins within the Active Region
            还有一个是峰值
            Peak fragment density of the Active Region (Value at “Summit”)
            我想要寻找在两组样本里差异最大的Active region,请问该怎么做呢?

            这是问题一 我的问题二是
            我还有一个表
            是所有的基因,还有它们在6个样本中的所有包含该gene的Active Region的峰值均值, 格式如

            NOC2L 12.25 13.75 13.75 14.75 13.50 18.00

            请问 想查看两组样本之间的差异基因,是不是可以用edgeR或者Deseq2分析了?如果可以,这个表达值是浮点数不是整数该怎么预处理一下?

          2. 问题1, 我倾向于使用均值。
            问题2, 我不建议使用edgeR或者Deseq。我建议使用anova test。

          3. 谢谢 按照您的说法 无论是我通过峰值均值找差异基因还是通过Average of fragment densities找出差异active region
            我都最好用anova做是吧?
            那么请问如何用anova做?
            是不是以下这样的思路:
            以这个基因为例
            NOC2L 条件1 12.25 13.75 13.75
            条件2 14.75 13.50 18.00
            进行两组的方差分析得到P值
            然后做一个循环 对对所有基因进行分析
            最后找出P<0.05的?

          4. 欧老师,
            能否请您看看我的代码看看对否?比如数据载入后还需要归一化么?e.g. voom(averagePeak,design,plot=TRUE,normalize=”quantile”)
            因为我算出来的adj.P有点怪,P值还可以。所以我不敢确定我的参数都是正确。

            我的代码:

            library(limma)
            averagePeak <- read.table("input.txt")
            SamplesCondition<-c(rep("C1",3),rep("C2",3))
            design <- model.matrix(~SamplesCondition)
            vmodel <- voom(averagePeak,design,plot=TRUE)
            fit <- lmFit(vmodel,design)
            fit <- eBayes(fit)
            topTable(fit,coef=2,sort="p")

            数据格式是这样的
            C1.R1 C1.R2 C1.R3 C2.R1 C2.R2 C2.R3
            Gene737 22.00 22.00 18.00 23.00 22.0 22.0
            Gene418 7.00 12.00 6.00 5.00 5.5 13.5
            Gene8069 85.00 114.00 94.00 123.00 129.0 106.0

            结果是
            logFC AveExpr t P.Value adj.P.Val B
            Gene8959 0.6314207 6.236774 7.429729 8.434147e-05 0.9124482 -2.760385
            Gene365 0.6314207 6.236774 7.429729 8.434147e-05 0.9124482 -2.760385

          5. 你的原始数据为啥会有小数?是normalize过的吗?如果是这样的话,voom是否合适?你应该仔细读一读voom的说明。还有,有些时候,adjust.p.value可能会看上去不好看,这只是一个BH重排的结果,所以必要时需要实验验证几个为准。

          6. 我在之前和您说过啊,小数是6个样本中的所有包含该gene的Active Region的峰值均值,没有normalize过。之前您推荐我不要用edgeR,用limma。
            adj.P确实不好看,都是0.9多,也是醉了。
            1)我的代码有无问题?
            2)您是建议我就用P值么?

          7. 我的建议是用limma,可是没让你用voom来做标准化啊。你直接把数据用library size修正后做log2转换后做lmFit就好了啊。

          8. 欧老师,如何按照您说的library size修正再做log2转换直接做lmfit?
            我看了limma的示例,我知道Lmfit(A,B)的B是一个model.matrix,但是A必须是一个data.object是吧?如何定义呢?我已经将数据用log2转换了,可是一载入Lmfit就报错啊。。。
            不知道您能不能把我的script改一改?谢谢。

  2. 老师您好,刚开始接触ChIP-Seq,想请问您一下,对于chip数据测序时,是选择单端测序比较好,还是双端测序较好,看到后面的peak detecter软件,大多都是支持单端测序,您在示例中提到的macs也是单端reads测序,查看下新版本的macs2已经支持了双端测序reads比对结果,想请教一下,两种测序方法对后续研究有什么影响。

    1. 我的观点是对于ChIP-Seq而言,测序深度是优先考虑的,单双端对于结果的影响小于测序深度的影响。在尽量经济的情况下,我会选择两倍深度的单端而非一倍深度的双端。双端测序对于定位可以map至多位点的reads有一定的好处,但是对于ChIP-Seq而言,好处有限。

  3. 老师您好,请教个问题,导师想做Chip-seq的组间差异分析,您有没有推荐的软件或者文章啊?

发表评论

电子邮件地址不会被公开。 必填项已用*标注