发布日期:2024-11-06 08:40 点击次数:137
图片
图片
图片
图片
图片
图片
今天学习下若何下载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.2150740MR分析
领先是加载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数据腹地处理):
半决赛前两局,第一盘柯洁获胜,占得先机;但在第2局比赛中执黑出战的柯洁由于中盘出现失误被一力辽扳回一局。
# 这一步需要联网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物联网app开发
本站仅提供存储处事,总计本体均由用户发布,如发现存害或侵权本体,请点击举报。