PAML之codeml正选择问题 2020-10-10

原理参考：http://www.chenlianfu.com/?p=3084
解析参考：https://user.qzone.qq.com/58001704?source=grouplist&t=0.07615071655538241
Branch models主要用于对系统发育树中不同支系 ω值差异性进行界定，主要有三个模型：

（1）One-ratio model (model = 0)：假设系统发育树中所有支系的 ω 值相等；

（2）Free-ratio model (model = 1 )：假设系统发育树中所有支系的 ω 值不相等；
注意此时的树未标定

outfile = branch.freeratio.mlc
     treefile = ../../tree/tree.raw

（3）Two-ratio model (model = 2)：假设前景枝和背景枝的ω 值不同；

seqtype = 1  * 1:codons; 2:AAs; 3:codons-->AAs
model = 1              * models for codons:
                       * 0:one, 1:b, 2:2 or more dN/dS ratios for branches

                       * models for AAs or codon-translated AAs:
                       * 0:poisson, 1:proportional, 2:Empirical, 3:Empirical+F
                       * 6:FromCodon, 7:AAClasses, 8:REVaa_0, 9:REVaa(nr=189)

image.png

零假设模型(null model)—one ratio model的 ω 值：omega (dN/dS) = 0.1
代表零假设下，整个树上的平均 ω 值为0.1
备择假设模型(branch model)—two ratio model的 ω 值：ω (dN/dS) for branches: 0.0999 0.22497
代表备择假设下，背景支的平均 ω 值为 0.0999 ，前景支的平均 ω 值为 0.22497。
若p值<0.01，接受备择假设，代表前景支的 ω 值明显比背景支大，由于前景支和背景支的 ω 值都<1。所以前景支应该经历了负选择的放松。

参考：https://yanzhongsino.github.io/2023/10/10/evolution_selection_paml.codeml_branch.model/

branch-site.fix.ctl

model = 2
                   * models for codons:
                       * 0:one, 1:b, 2:2 or more dN/dS ratios for branches
                   * models for AAs or codon-translated AAs:
                       * 0:poisson, 1:proportional, 2:Empirical, 3:Empirical+F
                       * 6:FromCodon, 7:AAClasses, 8:REVaa_0, 9:REVaa(nr=189)
      NSsites = 2    * 0:one w;1:neutral;2:selection; 3:discrete;4:freqs;
                   * 5:gamma;6:2gamma;7:beta;8:beta&w;9:beta&gamma;
                   * 10:beta&gamma+1; 11:beta&normal>1; 12:0&2normal>1;
                   * 13:3normal>0

        icode = 4  * 0:universal code; 1:mammalian mt; 2-10:see below
        Mgene = 0
                   * codon: 0:rates, 1:separate; 2:diff pi, 3:diff kapa, 4:all diff
                   * AA: 0:rates, 1:separate

    fix_kappa = 0  * 1: kappa fixed, 0: kappa to be estimated
        kappa = 2  * initial or fixed kappa
    fix_omega = 1    * 1: omega or omega_1 fixed, 0: estimate
        omega = 1   * initial or fixed omega, for codons or codon-based AAs

    fix_alpha = 1  * 0: estimate gamma shape parameter; 1: fix it at alpha
        alpha = 0. * initial or fixed alpha, 0:infinity (constant rate)
       Malpha = 0  * different alphas for genes
        ncatG = 8  * # of categories in dG of NSsites models

        getSE = 0  * 0: don't want them, 1: want S.E.s of estimates
 RateAncestor = 1  * (0,1,2): rates (alpha>0) or ancestral states (1 or 2)

   Small_Diff = .5e-6
    cleandata = 1  * remove sites with ambiguity data (1:yes, 0:no)?

branch-site.nofix.ctl

outfile = branch-site.DESI_ARAQ.nofix.mlc
     treefile = ../../tree/tree.DESI_ARAQ
    model = 2
                   * models for codons:
                       * 0:one, 1:b, 2:2 or more dN/dS ratios for branches
                   * models for AAs or codon-translated AAs:
                       * 0:poisson, 1:proportional, 2:Empirical, 3:Empirical+F
                       * 6:FromCodon, 7:AAClasses, 8:REVaa_0, 9:REVaa(nr=189)
      NSsites = 2    * 0:one w;1:neutral;2:selection; 3:discrete;4:freqs;
                   * 5:gamma;6:2gamma;7:beta;8:beta&w;9:beta&gamma;
                   * 10:beta&gamma+1; 11:beta&normal>1; 12:0&2normal>1;
                   * 13:3normal>0

        icode = 4  * 0:universal code; 1:mammalian mt; 2-10:see below
        Mgene = 0
                   * codon: 0:rates, 1:separate; 2:diff pi, 3:diff kapa, 4:all diff
                   * AA: 0:rates, 1:separate

    fix_kappa = 0  * 1: kappa fixed, 0: kappa to be estimated
        kappa = 2  * initial or fixed kappa
    fix_omega = 0     * 1: omega or omega_1 fixed, 0: estimate
        omega = 1.5    * initial or fixed omega, for codons or codon-based AAs

    fix_alpha = 1  * 0: estimate gamma shape parameter; 1: fix it at alpha
        alpha = 0. * initial or fixed alpha, 0:infinity (constant rate)
       Malpha = 0  * different alphas for genes
        ncatG = 8  * # of categories in dG of NSsites models

        getSE = 0  * 0: don't want them, 1: want S.E.s of estimates
 RateAncestor = 1  * (0,1,2): rates (alpha>0) or ancestral states (1 or 2)

   Small_Diff = .5e-6
    cleandata = 1  * remove sites with ambiguity data (1:yes, 0:no)?

最后编辑于：2024.05.23 12:34:45