免费视频淫片aa毛片_日韩高清在线亚洲专区vr_日韩大片免费观看视频播放_亚洲欧美国产精品完整版

打開APP
userphoto
未登錄

開通VIP,暢享免費電子書等14項超值服

開通VIP
一文搞定基因型數據清洗

大家好,我是飛哥。

我們知道很多分析之前,都要做基因型數據清洗,包括:

  • GWAS分析
  • GS分析
  • ……

這里介紹一下常用的基因型數據清洗方法。

  • 數據

  • 1 二進制文件

  • 2. plink二進制文件變?yōu)槲谋疚募╬ed和map)

  • 3. plink將vcf轉化為plink文件

  • 4. 提取樣本和SNP

    • 4.1 提取樣本

    • 4.2 提取SNP

  • 5. plink和表型數據合并

  • 6. 數據匯總

    • 6.1 次等位基因頻率(maf)

    • 6.2 缺失

    • 6.3 哈溫檢測

    • 6.4 雜合率檢測

  • 7. 數據質控qc

    • 7.1 樣本質控--mind

    • 7.2 位點質控--geno

    • 7.3 maf質控--maf

    • 7.4 哈溫質控 --hwe

數據

《統(tǒng)計遺傳學》中的章節(jié)介紹,有關代碼實操部分,單獨列出來,進行展示。

我已經下載整理好了,下載本書的電子版pdf+數據+代碼,鏈接:書籍及配套代碼領取--統(tǒng)計遺傳分析導論

1 二進制文件

文件中包括二進制的三個文件:

2. plink二進制文件變?yōu)槲谋疚募╬ed和map)

命令

 plink --bfile hapmap-ceu --recode --out hapmap-ceu
  • --bfile 是指定二進制plink的前綴名稱
  • --recode是生成文本ped和map的二進制文件
  • --out是指定輸出的結果文件前綴名稱

日志

PLINK v1.90b5.3 64-bit (21 Feb 2018)           www.cog-genomics.org/plink/1.9/
(C) 2005-2018 Shaun Purcell, Christopher Chang   GNU General Public License v3
Logging to hapmap-ceu.log.
Options in effect:
  --bfile hapmap-ceu
  --out hapmap-ceu
  --recode

1031523 MB RAM detected; reserving 515761 MB for main workspace.
2239392 variants loaded from .bim file.
60 people (30 males, 30 females) loaded from .fam.
Using 1 thread (no multithreaded calculations invoked).
Before main variant filters, 60 founders and 0 nonfounders present.
Calculating allele frequencies... done.
Warning: 59 het. haploid genotypes present (see hapmap-ceu.hh ); many commands
treat these as missing.
Total genotyping rate is 0.992022.
2239392 variants and 60 people pass filters and QC.
Note: No phenotypes present.
--recode ped to hapmap-ceu.ped + hapmap-ceu.map ... done.

上面的信息有:

  • 共有2239392個SNP位點
  • 共有60個個體,其中30個males,30個female
  • 結果生成hapmap-ceu.pedhapmap-ceu.map文件

3. plink將vcf轉化為plink文件

文件:ALL.chr21.vcf.gz

plink將vcf文件變?yōu)閜link的二進制文件(bed和bim和fam)。

命令:

plink --vcf ALL.chr21.vcf.gz --make-bed --out test_vcf

日志:

PLINK v1.90b5.3 64-bit (21 Feb 2018)           www.cog-genomics.org/plink/1.9/
(C) 2005-2018 Shaun Purcell, Christopher Chang   GNU General Public License v3
Logging to test_vcf.log.
Options in effect:
  --make-bed
  --out test_vcf
  --vcf ALL.chr21.vcf.gz

1031523 MB RAM detected; reserving 515761 MB for main workspace.
--vcf: test_vcf-temporary.bed + test_vcf-temporary.bim + test_vcf-temporary.fam
written.
1105538 variants loaded from .bim file.
2504 people (0 males, 0 females, 2504 ambiguous) loaded from .fam.
Ambiguous sex IDs written to test_vcf.nosex .
Using 1 thread (no multithreaded calculations invoked).
Before main variant filters, 2504 founders and 0 nonfounders present.
Calculating allele frequencies... done.
Total genotyping rate is 0.999787.
1105538 variants and 2504 people pass filters and QC.
Note: No phenotypes present.
--make-bed to test_vcf.bed + test_vcf.bim + test_vcf.fam ... done.

4. 提取樣本和SNP

  • 可以使用--keep 和 --remove 進行樣本ID的提取或者刪除
  • 可以使用--extract和 --exclude進行SNP的提取或者刪除

樣本ID的示例:兩列FID和IID,沒有行頭:

SNP的ID的示例:一列SNP名稱,沒有行頭:

4.1 提取樣本

代碼:

 plink --bfile hapmap-ceu --keep list.txt --make-bed --out selectedIndividuals

日志:

PLINK v1.90b5.3 64-bit (21 Feb 2018)           www.cog-genomics.org/plink/1.9/
(C) 2005-2018 Shaun Purcell, Christopher Chang   GNU General Public License v3
Logging to selectedIndividuals.log.
Options in effect:
  --bfile hapmap-ceu
  --keep list.txt
  --make-bed
  --out selectedIndividuals

1031523 MB RAM detected; reserving 515761 MB for main workspace.
2239392 variants loaded from .bim file.
60 people (30 males, 30 females) loaded from .fam.
--keep: 60 people remaining.
Using 1 thread (no multithreaded calculations invoked).
Before main variant filters, 60 founders and 0 nonfounders present.
Calculating allele frequencies... done.
Warning: 59 het. haploid genotypes present (see selectedIndividuals.hh ); many
commands treat these as missing.
Total genotyping rate is 0.992022.
2239392 variants and 60 people pass filters and QC.
Note: No phenotypes present.
--make-bed to selectedIndividuals.bed + selectedIndividuals.bim +
selectedIndividuals.fam ... done.

4.2 提取SNP

代碼:

plink --bfile hapmap-ceu --extract snp_name.txt --make-bed --out selectedSNP

日志:

PLINK v1.90b5.3 64-bit (21 Feb 2018)           www.cog-genomics.org/plink/1.9/
(C) 2005-2018 Shaun Purcell, Christopher Chang   GNU General Public License v3
Logging to selectedSNP.log.
Options in effect:
  --bfile hapmap-ceu
  --extract snp_name.txt
  --make-bed
  --out selectedSNP

1031523 MB RAM detected; reserving 515761 MB for main workspace.
2239392 variants loaded from .bim file.
60 people (30 males, 30 females) loaded from .fam.
--extract: 20 variants remaining.
Using 1 thread (no multithreaded calculations invoked).
Before main variant filters, 60 founders and 0 nonfounders present.
Calculating allele frequencies... done.
Total genotyping rate is 0.995833.
20 variants and 60 people pass filters and QC.
Note: No phenotypes present.
--make-bed to selectedSNP.bed + selectedSNP.bim + selectedSNP.fam ... done.

也可以提取單個SNP:

 plink --file hapmap-ceu --snps rs9930506 --make-bed --out rs9930506sample

結果:

 rs9930506sample.bed + rs9930506sample.bim + rs9930506sample.fam

5. plink和表型數據合并

如果想要把表型數據和基因型數據合并,需要整理的表型格式:FID,IID,y三列。

FID IID BMI
0 HG00096 25.022827
0 HG00097 24.853638
0 HG00099 23.689295
0 HG00100 27.016203
0 HG00101 21.461624
0 HG00102 20.673635
0 HG00103 25.71508
0 HG00104 25.252243
0 HG00106 22.765049

看一下原始的fam文件:

0 HG00096 0 0 1 -9
0 HG00097 0 0 2 -9
0 HG00099 0 0 2 -9
0 HG00100 0 0 2 -9
0 HG00101 0 0 1 -9
0 HG00102 0 0 2 -9
0 HG00103 0 0 1 -9
0 HG00104 0 0 2 -9
0 HG00106 0 0 2 -9
0 HG00108 0 0 1 -9

合并命令:

 plink --bfile 1kg_EU_qc --pheno BMI_pheno.txt --make-bed --out 1kg_EU_BMI

生成文件:

1kg_EU_BMI.bed + 1kg_EU_BMI.bim + 1kg_EU_BMI.fam

看一下新的fam文件:

0 HG00096 0 0 1 25.0228
0 HG00097 0 0 2 24.8536
0 HG00099 0 0 2 23.6893
0 HG00100 0 0 2 27.0162
0 HG00101 0 0 1 21.4616
0 HG00102 0 0 2 20.6736
0 HG00103 0 0 1 25.7151
0 HG00104 0 0 2 25.2522
0 HG00106 0 0 2 22.765
0 HG00108 0 0 1 23.069

可以看到,已經合并成功了。

6. 數據匯總

6.1 次等位基因頻率(maf)

查看基因頻率的統(tǒng)計結果,用--freq

命令:

plink --bfile hapmap-ceu --freq --out Allele_Frequency

日志:

PLINK v1.90b6.21 64-bit (19 Oct 2020)          www.cog-genomics.org/plink/1.9/
(C) 2005-2020 Shaun Purcell, Christopher Chang   GNU General Public License v3
Logging to Allele_Frequency.log.
Options in effect:
  --bfile hapmap-ceu
  --freq
  --out Allele_Frequency

15236 MB RAM detected; reserving 7618 MB for main workspace.
2239392 variants loaded from .bim file.
60 people (30 males, 30 females) loaded from .fam.
Using 1 thread (no multithreaded calculations invoked).
Before main variant filters, 60 founders and 0 nonfounders present.
Calculating allele frequencies... done.
Total genotyping rate is 0.992022.
--freq: Allele frequencies (founders only) written to Allele_Frequency.frq .
Warning: 59 het. haploid genotypes present (see Allele_Frequency.hh )
; many
commands treat these as missing.

查看結果,結果文件是Allele_Frequency.frq

$ head Allele_Frequency.frq
 CHR          SNP   A1   A2          MAF  NCHROBS
   1   rs12565286    C    G       0.0678      118
   1   rs12138618    A    G      0.05833      120
   1    rs3094315    G    A       0.1552      116
   1    rs3131968    A    G        0.125      120
   1   rs12562034    A    G      0.09167      120
   1    rs2905035    A    G        0.125      120
   1   rs12124819    G    A       0.3417      120
   1    rs2980319    A    T        0.125      120
   1    rs4040617    G    A        0.125      120

如果對其進行質控,用--maf 0.01,會去掉maf小于0.01的位點。

6.2 缺失

缺失包括樣本缺失率統(tǒng)計和位點缺失率統(tǒng)計。

命令:

plink --bfile hapmap-ceu --missing --out missing_data

日志:

$ plink --bfile hapmap-ceu --missing --out missing_data
PLINK v1.90b6.21 64-bit (19 Oct 2020)          www.cog-genomics.org/plink/1.9/
(C) 2005-2020 Shaun Purcell, Christopher Chang   GNU General Public License v3
Logging to missing_data.log.
Options in effect:
  --bfile hapmap-ceu
  --missing
  --out missing_data

15236 MB RAM detected; reserving 7618 MB for main workspace.
2239392 variants loaded from .bim file.
60 people (30 males, 30 females) loaded from .fam.
Using 1 thread (no multithreaded calculations invoked).
Before main variant filters, 60 founders and 0 nonfounders present.
Calculating allele frequencies... done.
Total genotyping rate is 0.992022.
--missing: Sample missing data report written to missing_data.imiss, and
variant-based missing data report written to missing_data.lmiss.
Warning: 59 het. haploid genotypes present (see missing_data.hh ); many
commands treat these as missing.

查看樣本缺失的文件:

$ head missing_data.imiss
 FID       IID MISS_PHENO   N_MISS   N_GENO   F_MISS
1334   NA12144          Y    15077  2239392 0.006733
1334   NA12145          Y    19791  2239392 0.008838
1334   NA12146          Y    13981  2239392 0.006243
1334   NA12239          Y    14072  2239392 0.006284
1340   NA06994          Y    16080  2239392 0.007181
1340   NA07000          Y    26113  2239392  0.01166
1340   NA07022          Y    17467  2239392   0.0078
1340   NA07056          Y    12133  2239392 0.005418
1341   NA07034          Y    20425  2239392 0.009121

查看位點缺失的文件:

$ head missing_data.lmiss
 CHR          SNP   N_MISS   N_GENO   F_MISS
   1   rs12565286        1       60  0.01667
   1   rs12138618        0       60        0
   1    rs3094315        2       60  0.03333
   1    rs3131968        0       60        0
   1   rs12562034        0       60        0
   1    rs2905035        0       60        0
   1   rs12124819        0       60        0
   1    rs2980319        0       60        0
   1    rs4040617        0       60        0

6.3 哈溫檢測

--hardy

代碼:

plink --bfile 1kg_EU_qc --hardy

日志:

PLINK v1.90b5.3 64-bit (21 Feb 2018)           www.cog-genomics.org/plink/1.9/
(C) 2005-2018 Shaun Purcell, Christopher Chang   GNU General Public License v3
Logging to plink.log.
Options in effect:
  --bfile 1kg_EU_qc
  --hardy

1031523 MB RAM detected; reserving 515761 MB for main workspace.
851065 variants loaded from .bim file.
379 people (178 males, 201 females) loaded from .fam.
Using 1 thread (no multithreaded calculations invoked).
Before main variant filters, 379 founders and 0 nonfounders present.
Calculating allele frequencies... done.
--hardy: Writing Hardy-Weinberg report (founders only) to plink.hwe ... done.

結果:

$ head plink.hwe
 CHR         SNP     TEST   A1   A2                 GENO   O(HET)   E(HET)            P
   1   rs1048488  ALL(NP)    C    T            9/136/234   0.3588   0.3238      0.03894
   1   rs3115850  ALL(NP)    T    C            8/135/236   0.3562    0.319      0.02412
   1   rs2519031  ALL(NP)    G    A             0/11/368  0.02902   0.0286            1
   1   rs4970383  ALL(NP)    A    C           22/149/208   0.3931   0.3796       0.5881
   1   rs4475691  ALL(NP)    T    C           13/118/248   0.3113   0.3078            1
   1   rs1806509  ALL(NP)    C    A           67/171/141   0.4512   0.4809       0.2402
   1   rs7537756  ALL(NP)    G    A           14/124/241   0.3272   0.3206       0.8725
   1  rs28576697  ALL(NP)    C    T           37/161/181   0.4248   0.4278       0.9045
   1   rs7523549  ALL(NP)    T    C             0/28/351  0.07388  0.07115            1
[dengfei@localhost 08_chapter]$ wc -l plink.hwe
851066 plink.hwe

6.4 雜合率檢測

--het

代碼:

plink --bfile 1kg_EU_qc --het

日志:

PLINK v1.90b5.3 64-bit (21 Feb 2018)           www.cog-genomics.org/plink/1.9/
(C) 2005-2018 Shaun Purcell, Christopher Chang   GNU General Public License v3
Logging to plink.log.
Options in effect:
  --bfile 1kg_EU_qc
  --het

1031523 MB RAM detected; reserving 515761 MB for main workspace.
851065 variants loaded from .bim file.
379 people (178 males, 201 females) loaded from .fam.
Using 1 thread (no multithreaded calculations invoked).
Before main variant filters, 379 founders and 0 nonfounders present.
Calculating allele frequencies... done.
851065 variants and 379 people pass filters and QC.
Note: No phenotypes present.
--het: 851065 variants scanned, report written to plink.het .


結果:

$ head plink.het
 FID       IID       O(HOM)       E(HOM)        N(NM)            F
   0   HG00096       582990    5.838e+05       851065    -0.003127
   0   HG00097       582246    5.838e+05       851065    -0.005911
   0   HG00099       582456    5.838e+05       851065    -0.005125
   0   HG00100       582700    5.838e+05       851065    -0.004212
   0   HG00101       581527    5.838e+05       851065    -0.008602
   0   HG00102       585010    5.838e+05       851065     0.004432
   0   HG00103       583984    5.838e+05       851065    0.0005925
   0   HG00104       586189    5.838e+05       851065     0.008844
   0   HG00106       583868    5.838e+05       851065    0.0001584
$ wc -l plink.het
380 plink.het


7. 數據質控qc

7.1 樣本質控--mind

代碼:

 plink --bfile 1kg_EU_qc --mind 0.1 --recode --out re1

日志:


PLINK v1.90b5.3 64-bit (21 Feb 2018)           www.cog-genomics.org/plink/1.9/
(C) 2005-2018 Shaun Purcell, Christopher Chang   GNU General Public License v3
Logging to re1.log.
Options in effect:
  --bfile 1kg_EU_qc
  --mind 0.1
  --out re1
  --recode

1031523 MB RAM detected; reserving 515761 MB for main workspace.
851065 variants loaded from .bim file.
379 people (178 males, 201 females) loaded from .fam.
0 people removed due to missing genotype data (--mind).
Using 1 thread (no multithreaded calculations invoked).
Before main variant filters, 379 founders and 0 nonfounders present.
Calculating allele frequencies... done.
851065 variants and 379 people pass filters and QC.
Note: No phenotypes present.
--recode ped to re1.ped + re1.map ... done.


結果:


$ wc -l re1.*
        24 re1.log
    851065 re1.map
       379 re1.ped
    851468 總用量

7.2 位點質控--geno

代碼:

plink --bfile 1kg_EU_qc --geno 0.1 --recode --out re1

日志:

PLINK v1.90b5.3 64-bit (21 Feb 2018)           www.cog-genomics.org/plink/1.9/
(C) 2005-2018 Shaun Purcell, Christopher Chang   GNU General Public License v3
Logging to re1.log.
Options in effect:
  --bfile 1kg_EU_qc
  --geno 0.1
  --out re1
  --recode

1031523 MB RAM detected; reserving 515761 MB for main workspace.
851065 variants loaded from .bim file.
379 people (178 males, 201 females) loaded from .fam.
Using 1 thread (no multithreaded calculations invoked).
Before main variant filters, 379 founders and 0 nonfounders present.
Calculating allele frequencies... done.
0 variants removed due to missing genotype data (--geno).
851065 variants and 379 people pass filters and QC.
Note: No phenotypes present.
--recode ped to re1.ped + re1.map ... done.



結果:

$ wc -l re1*
        24 re1.log
    851065 re1.map
       379 re1.ped
    851468 總用量

7.3 maf質控--maf

代碼:

plink --bfile 1kg_EU_qc --maf 0.01 --recode --out re1

日志:

PLINK v1.90b5.3 64-bit (21 Feb 2018)           www.cog-genomics.org/plink/1.9/
(C) 2005-2018 Shaun Purcell, Christopher Chang   GNU General Public License v3
Logging to re1.log.
Options in effect:
  --bfile 1kg_EU_qc
  --maf 0.01
  --out re1
  --recode

1031523 MB RAM detected; reserving 515761 MB for main workspace.
851065 variants loaded from .bim file.
379 people (178 males, 201 females) loaded from .fam.
Using 1 thread (no multithreaded calculations invoked).
Before main variant filters, 379 founders and 0 nonfounders present.
Calculating allele frequencies... done.
17767 variants removed due to minor allele threshold(s)
(--maf/--max-maf/--mac/--max-mac).
833298 variants and 379 people pass filters and QC.
Note: No phenotypes present.
--recode ped to re1.ped + re1.map ... done.



結果:

$ wc -l re1*
        25 re1.log
    833298 re1.map
       379 re1.ped
    833702 總用量


7.4 哈溫質控 --hwe

代碼:

plink --bfile 1kg_EU_qc --hwe 1e-5 --recode --out re1

日志:

PLINK v1.90b5.3 64-bit (21 Feb 2018)           www.cog-genomics.org/plink/1.9/
(C) 2005-2018 Shaun Purcell, Christopher Chang   GNU General Public License v3
Logging to re1.log.
Options in effect:
  --bfile 1kg_EU_qc
  --hwe 1e-5
  --out re1
  --recode

1031523 MB RAM detected; reserving 515761 MB for main workspace.
851065 variants loaded from .bim file.
379 people (178 males, 201 females) loaded from .fam.
Using 1 thread (no multithreaded calculations invoked).
Before main variant filters, 379 founders and 0 nonfounders present.
Calculating allele frequencies... done.
--hwe: 25 variants removed due to Hardy-Weinberg exact test.
851040 variants and 379 people pass filters and QC.
Note: No phenotypes present.
--recode ped to re1.ped + re1.map ... done.



結果:

$ wc -l re1*
        24 re1.log
    851040 re1.map
       379 re1.ped
    851443 總用量
本站僅提供存儲服務,所有內容均由用戶發(fā)布,如發(fā)現(xiàn)有害或侵權內容,請點擊舉報
打開APP,閱讀全文并永久保存 查看更多類似文章
猜你喜歡
類似文章
筆記GWAS 操作流程2-3:MAF過濾
生信小工具:Plink之常用操作實戰(zhàn)演練(2)
plink格式文件的介紹及相關轉換
plink合并文件并更新SNP位置(merge file, update SNP position)...
GWAS | 原理和流程 | 全基因組關聯(lián)分析 | Linkage disequilibrium (LD)連鎖不平衡 | 曼哈頓圖 Manhattan
GWAS Catalog數據庫簡介
更多類似文章 >>
生活服務
分享 收藏 導長圖 關注 下載文章
綁定賬號成功
后續(xù)可登錄賬號暢享VIP特權!
如果VIP功能使用有故障,
可點擊這里聯(lián)系客服!

聯(lián)系客服