图片

图片

图片

图片

图片

图片

今天学习下如何下载UKBiobank的数据并进行整理，把它酿成TwoSampleMR草率使用的法子。

数据下载

掀开ukbiobank的官网：https://www.nealelab.is/uk-biobank

掀开后往下拉，找到results can be found here，点击它：

图片

点击之后即可来到数据界面，率先是README这个sheet，这里会有一些先容，每个文献包含哪些实际，每个文献的列名是什么意旨道理等，皆有详备的讲解，各人一定要看！

图片

图片

为了便捷以后使用，这个文献是不错免费下载的，提倡各人径直下来，背面就不错用EXCEL随时掀开稽查了。

图片

点击底下的Manifest 201807就不错来到数据细目和下载界面了，UKbiobank这个免费的数据照旧好久不更新了，不错看到数据就留在2018年了。

图片

数据细目界面有每个数据的phenotype description、phenotype code、数据下载相连等信息。

不错看到每个表型（phenotype）皆有6个数据，分袂是inrt和raw各有3种（分为female、male、both_sexes）。

图片

把底下的动荡条往右边拉到临了，就不错看到每个数据的下载相连了，你鼠标放上去就会有复制相连的辅导，径直复制在浏览器或者迅雷等软件中掀开即可下载数据了。

每个表型皆有6个数据，那么咱们应该是用哪个呢？官方的github[1]先容说选拔irnt手脚接下来的分析：

Overall， the results are largely consistent regardless of the choice of IRNT(inverse rank-normal transformed) or raw untransformed phenotypes. Since IRNT does appear to provide a marginal boost to the h2g results， especially in terms of significance， we choose to treat the IRNT version as the primary analysis for continuous phenotypes. (Results for the raw， untransformed versions will still appear in the complete results file though.)

是以咱们也不纠结了，径直用irnt的数据分析，要是你的数据需要分开性别，那么你就单独下载某个性别的数据，要么即是下载both_sexes的。

数据读取及整理

我这里选拔了fat这个表型的数据（上头有展示）。下载好之后径直用vroom读取即可：

fat_bothsex <- vroom::vroom("./ukbiobank/100004_irnt.gwas.imputed_v3.both_sexes.tsv.bgz")dim(fat_bothsex)## [1] 13791467       11colnames(fat_bothsex)##  [1] "variant"                "minor_allele"           "minor_AF"              ##  [4] "low_confidence_variant" "n_complete_samples"     "AC"                    ##  [7] "ytx"                    "beta"                   "se"                    ## [10] "tstat"                  "pval"

这个数据亦然很大的，1000多万行，11列。

# 看下前6行head(fat_bothsex)## # A tibble: 6 × 11##   variant minor_allele minor_AF low_confidence_variant n_complete_samples     AC##   <chr>   <chr>           <dbl> <lgl>                               <dbl>  <dbl>## 1 1:1579… T             0       TRUE                                51453 0     ## 2 1:6948… A             2.04e-5 TRUE                                51453 2.10e0## 3 1:6956… C             1.80e-4 TRUE                                51453 1.85e1## 4 1:1398… T             2.03e-5 TRUE                                51453 2.09e0## 5 1:6927… C             1.12e-1 FALSE                               51453 1.15e4## 6 1:6937… G             1.17e-1 FALSE                               51453 1.20e4## # ℹ 5 more variables: ytx <dbl>， beta <dbl>， se <dbl>， tstat <dbl>， pval <dbl>

这个数据的每一列是什么意旨道理，皆在前边说过的README部分有详备的先容，各人一定要去看！

比如第一列是variant，它的意旨道理是chr:pos:ref:alt，alt是effect allele，需要使用这个和疑望文献归并。这些信息皆在README内部：

图片

这个数据是莫得rsid的，然则你仔细读过README之后就会发现，在variants.tsv.bgz这个文献中是有rsid的，况兼这个文献也有chr:pos:ref:alt这些信息：

图片

情欲印象下载

是以把这个文献下载下来，和咱们的fat_bothsex文献归并一下，就有rsid了。variants.tsv.bgz这个文献就在Manifest 201807的第11行，找到下载相连径直下载即可：

图片

下载好之后咱们把它读取进来：

variants_anno <- vroom::vroom("./ukbiobank/variants.tsv.bgz")dim(variants_anno)## [1] 13791467       25

发现这个数据和fat_bothsex这个数据的行数是同样的。其着实README的描摹中东说念主家就说的很线路了，这个律例即是同样的，你不错通过variant这一列进行归并，也不错径直复制粘贴到沿途：

图片

稽查前6行：

variants_anno[1:6，]## # A tibble: 6 × 25##   variant  chr      pos ref   alt   rsid  varid consequence consequence_category##   <chr>    <chr>  <dbl> <chr> <chr> <chr> <chr> <chr>       <chr>               ## 1 1:15791… 1      15791 C     T     rs54… 1:15… splice_reg… missense            ## 2 1:69487… 1      69487 G     A     rs56… 1:69… missense_v… missense            ## 3 1:69569… 1      69569 T     C     rs25… 1:69… missense_v… missense            ## 4 1:13985… 1     139853 C     T     rs53… 1:13… splice_reg… missense            ## 5 1:69279… 1     692794 CA    C     1:69… 1:69… upstream_g… non_coding          ## 6 1:69373… 1     693731 A     G     rs12… 1:69… upstream_g… non_coding          ## # ℹ 16 more variables: info <dbl>， call_rate <dbl>， AC <dbl>， AF <dbl>，## #   minor_allele <chr>， minor_AF <dbl>， p_hwe <dbl>， n_called <dbl>，## #   n_not_called <dbl>， n_hom_ref <dbl>， n_het <dbl>， n_hom_var <dbl>，## #   n_non_ref <dbl>， r_heterozygosity <dbl>， r_het_hom_var <dbl>，## #   r_expected_het_frequency <dbl>

望望是不是十足同样：

# 战胜是十足一致identical(variants_anno$variant， fat_bothsex$variant)## [1] TRUE

是以咱们径直归并一下即可：

fat_bothsex <- cbind.data.frame(variants_anno[，2:6]， fat_bothsex)head(fat_bothsex)##   chr    pos ref alt          rsid       variant minor_allele    minor_AF## 1   1  15791   C   T   rs547522712   1:15791:C:T            T 0.00000e+00## 2   1  69487   G   A   rs568226429   1:69487:G:A            A 2.03879e-05## 3   1  69569   T   C     rs2531267   1:69569:T:C            C 1.80062e-04## 4   1 139853   C   T   rs533633326  1:139853:C:T            T 2.03117e-05## 5   1 692794  CA   C 1:692794_CA_C 1:692794:CA:C            C 1.12102e-01## 6   1 693731   A   G    rs12238997  1:693731:A:G            G 1.16788e-01##   low_confidence_variant n_complete_samples          AC        ytx       beta## 1                   TRUE              51453     0.00000    0.00000        NaN## 2                   TRUE              51453     2.09804   -1.99658 -0.9567220## 3                   TRUE              51453    18.52940   -2.86950 -0.1222020## 4                   TRUE              51453     2.09020   -2.00107 -0.9583010## 5                  FALSE              51453 11535.90000 -162.96100 -0.0233942## 6                  FALSE              51453 12018.20000 -105.14000 -0.0126221##          se     tstat      pval## 1       NaN       NaN       NaN## 2 0.6939280 -1.378710 0.1679920## 3 0.2382830 -0.512845 0.6080620## 4 0.6939330 -1.380970 0.1672940## 5 0.0107131 -2.183700 0.0289888## 6 0.0101812 -1.239750 0.2150740

MR分析

率先是加载R包：

library(TwoSampleMR)

当前fat_bothsex这个数据是一个GWAS的原始数据，你不错把它用作是结局有关的数据，也不错把它用作是清楚有关的数据。

要是你要把这个数据用作是结局有关的数据，那即是径直进行以下代码即可：

# 扎眼列名逐一双应outcome_data <- format_data(dat = fat_bothsex，                      type = "outcome"， # 类型是outcome                      snp_col = "rsid"，                      beta_col = "beta"，                      se_col = "se"，                      eaf_col = "minor_AF"，                      effect_allele_col = "alt"，                      other_allele_col = "ref"，                      pval_col = "pval"，                      samplesize_col = "n_complete_samples"，                      #ncase_col = "ncase_col"，                      #ncontrol_col = "ncontrol_col"，                      #gene_col = "nearest_genes"，                      chr_col = "chr"，                      pos_col = "pos"                      )

当前这个outcome_data即是一个结局数据了，要是你照旧有了清楚数据，那么就不错径直进行MR分析了。

要是你要把这个数据用作是清楚数据，那么率先需要证实P值筛选一下：

# 先去掉NaNdf1 <- fat_bothsex[!is.nan(fat_bothsex$pval)，]# 证实P值筛选，我这个数据莫得小于5e-8的，是以选拔了5e-6手脚演示df1 <- subset(df1， pval < 5e-6)

然后再法子化数据：

# 扎眼列名逐一双应exposure_data <- format_data(dat = df1，                      type = "exposure"， # 类型是exposure                      snp_col = "rsid"，                      beta_col = "beta"，                      se_col = "se"，                      eaf_col = "minor_AF"，                      effect_allele_col = "alt"，                      other_allele_col = "ref"，                      pval_col = "pval"，                      samplesize_col = "n_complete_samples"，                      #ncase_col = "ncase_col"，                      #ncontrol_col = "ncontrol_col"，                      #gene_col = "nearest_genes"，                      chr_col = "chr"，                      pos_col = "pos"                      )

然后再开动下clump（这个经由需要联网的，要是你思十足土产货化开动，可参考下载IEU OPEN GWAS数据土产货处理）：

# 这一步需要联网exposure_data <- clump_data(exposure_data，                            clump_kb=10000， clump_r2=0.001)

当前这个exposure_data即是一个清楚数据了，要是你照旧有了结局数据，那么就不错径直进行MR分析了。

其中还有一些问题并莫得说，比如一个ref有多个alt，或者minor_allele有多个等问题，以后再先容，各人遭受了也不错我方搜索下。

贬责！easy！

参考贵寓[1]

github: https://nealelab.github.io/UKBB_ldsc/select_topline.html文爱 电报群

本站仅提供存储就业，整个实际均由用户发布，如发现存害或侵权实际，请点击举报。