乡下人产国偷v产偷v自拍,国产午夜片在线观看,婷婷成人亚洲综合国产麻豆,久久综合给合久久狠狠狠9

  • <output id="e9wm2"></output>
    <s id="e9wm2"><nobr id="e9wm2"><ins id="e9wm2"></ins></nobr></s>

    • 分享

      筆記 | GWAS 操作流程5-1:根紅苗正的GWAS分析軟件:GEMMA

       育種數(shù)據(jù)分析 2021-11-18

      筆記 GWAS 操作流程5-1:根紅苗正的GWAS分析軟件:GEMMA

      1. GEMMA軟件介紹

      這個(gè)肯定厲害了,是「大家閨秀」,是「名門(mén)望族」,是「根紅苗正」的GWAS分析軟件。

      「GEMMA名稱(chēng)來(lái)源:」

      • G:Genome-wide
      • E:Efficient
      • MM:Mixed-model
      • A:Association

      「GEMMAX主要特點(diǎn):」

      ?

      快,話(huà)說(shuō)同樣的檢測(cè)方法,GEMMA跑了3.3小時(shí),而EMMA估計(jì)要跑27年???

      ?

      2. GEMMA語(yǔ)法特點(diǎn)

      相對(duì)于plink的語(yǔ)法,GEMMA語(yǔ)法更簡(jiǎn)練,一個(gè)杠,一個(gè)字母。比如:

      • 表型數(shù)據(jù):-p
      • 協(xié)變量:-c而plink的語(yǔ)法的是兩個(gè)杠,一個(gè)單詞:
      • 表型數(shù)據(jù):--pheno
      • 協(xié)變量:--covar

      GEMMA支持plink的二進(jìn)制文件:

      • 讀取plink文件:-bfile

      「表型數(shù)據(jù)格式:」一列,注意順序和基因組的ID順序一致,如果是多個(gè)性狀,那就是多列,沒(méi)有ID列。

      「協(xié)變量格式:」協(xié)變量是數(shù)字類(lèi)型,如果有因子,需要轉(zhuǎn)化為虛擬變量的形式。沒(méi)有ID列,第一列是cov1,第二列是cov2……,注意,如果有一協(xié)變量,第一列為1(截距),第二列為協(xié)變量。

      比如下面的數(shù)據(jù)是一個(gè)協(xié)變量,第一列為截距。

      3. GEMMA分析一般線(xiàn)性模型沒(méi)有協(xié)變量

      「首先將plink格式轉(zhuǎn)化為二進(jìn)制的plink格式:」

      plink --file b --make-bed --out c

      「然后將表型數(shù)據(jù)提取單獨(dú)一列:」

      awk '{print $3}' phe.txt >p.txt

      「然后進(jìn)行一般線(xiàn)性模型關(guān)聯(lián)分析:」

      gemma-0.98.1-linux-static -bfile c -p p.txt -lm 1

      「結(jié)果和plink的linear結(jié)果對(duì)比:」

      plink的結(jié)果:gemma的結(jié)果:

      兩者結(jié)果完全一致。

      事實(shí)上,加上協(xié)變量的分析,gemma和plink的結(jié)果也是一樣的,因?yàn)槎际菓?yīng)用的是一般線(xiàn)性模型。

      4. GEMMA分析混合線(xiàn)性模型

      「第一步:先生成G矩陣」

      gemma-0.98.1-linux-static -bfile c -gk 2 -p p.txt 

      代碼解釋?zhuān)?/p>

      • -bfile 讀取plink的二進(jìn)制文件
      • -gk 2 標(biāo)準(zhǔn)化的方法計(jì)算G矩陣
      • -p 讀取表型數(shù)據(jù)(這個(gè)不能省掉)
      GEMMA 0.98.1 (2018-12-10by Xiang Zhou and team (C) 2012-2018
      Reading Files ... 
      ## number of total individuals = 1500
      ## number of analyzed individuals = 1500
      ## number of covariates = 1
      ## number of phenotypes = 1
      ## number of total SNPs/var        =    10000
      ## number of analyzed SNPs         =     3946
      Calculating Relatedness Matrix ... 
      ================================================== 100%
      **** INFO: Done.

      G矩陣在output文件夾下:result.sXX.txt

      「第二步:使用混合線(xiàn)性模型進(jìn)行GWAS分析」

      gemma-0.98.1-linux-static -bfile c -k output/result.sXX.txt -lmm 1 -p p.txt 
      代碼解釋?zhuān)?br>* -bfile plink的二進(jìn)制文件
      * -k 讀取G矩陣的文件
      * -lmm 1 使用Wald的方法進(jìn)行SNP檢驗(yàn)
      * -p 表型數(shù)據(jù)

      GEMMA 0.98.1 (2018-12-10by Xiang Zhou and team (C) 2012-2018
      Reading Files ... 
      ## number of total individuals = 1500
      ## number of analyzed individuals = 1500
      ## number of covariates = 1
      ## number of phenotypes = 1
      ## number of total SNPs/var        =    10000
      ## number of analyzed SNPs         =     3946
      Start Eigen-Decomposition...
      pve estimate 
      =0.120599
      se(pve) =0.0279188
      ================================================== 100%
      **** INFO: Done.

      「第三步:查看結(jié)果文件」結(jié)果在output文件夾下:result.assoc.txt

      5. GEMMA中LM模型和LMM模型的結(jié)果比較

      setwd("/home/dengfei/gwas/qmsim/dat/plink_file/10_gemma_analysis_lmm/output")

      mm_re = fread("result.assoc.txt")
      head(mm_re)

      lm_re = fread("../../06_gemma_analysis_lm/output/result.assoc.txt")
      head(lm_re)
      lm_re1 = lm_re[!is.na(p_wald)]
      head(lm_re1)

      dim(mm_re)
      dim(lm_re1)

      head(mm_re)
      head(lm_re1)

      re1 = merge(mm_re,lm_re1,by="rs")
      head(re1)

      # Pvalue 比較
      cor(re1$p_wald.x,re1$p_wald.y)
      plot(re1$p_wald.x,re1$p_wald.y)

      # Beta回歸系數(shù)比較
      cor(re1$beta.x,re1$beta.y)
      plot(re1$beta.x,re1$beta.y)

      「Pvalue比較」

      > cor(re1$p_wald.x,re1$p_wald.y)
      [10.5278895

      「Beta回歸系數(shù)比較:」

      > cor(re1$beta.x,re1$beta.y)
      [10.8357564

      「注意:」這里面,混合線(xiàn)性模型分析時(shí),一般要加上PCA的結(jié)果,防止群體結(jié)構(gòu)造成的假陽(yáng)性。下一章節(jié),我們介紹LMM模型,如何加入?yún)f(xié)變量。

      寫(xiě)到這里,流暢的感覺(jué)撲鼻而來(lái),不寫(xiě)了,發(fā)個(gè)朋友圈,洗洗睡吧……

        轉(zhuǎn)藏 分享 獻(xiàn)花(0

        0條評(píng)論

        發(fā)表

        請(qǐng)遵守用戶(hù) 評(píng)論公約

        類(lèi)似文章 更多