前言

开始本节之前，可先预先查看标准数据格式，当然也可以自己按照标准格式准备输入文件。

在开始分析之前，都需要进行数据预处理，将所有数据都准备成统一数据格式，方便后续分析。因此从原始数据的标准化格式准备入手！

如果用常规流程处理的话，往往需要非常大的计算机内存，并且耗时很长。而且文件格式不尽相同，所以为了满足内存占用少，耗时短；研究了常用的几大数据库文件内容，并开发了内置的文件转换工具，已内置到R包。

所以接下来，我们按照数据来源进行分类，对GWAS原始数据进行标准化准备！

IEU数据
GWAS Catalog
UKB数据
FinnGen
其他来源数据
deCODE2021 血浆蛋白
UKB-PPP2023 血浆蛋白

IEU数据

示例处理ieu-a-302.vcf.gz，源文件的大小是61.3MB。

# 加载包
library(DrugTargetMR)

# 调用函数
DrugTargetMR::prepare_ieu(file_path = "ieu-a-1126.vcf.gz",
                          out_path="./",
                          generate_mr = T,
                          generate_smr = T)

# 2023/11/17 15:13:12 IEU数据文件路径： ./ieu-a-302.vcf.gz
# 2023/11/17 15:13:12 输出mr文件路径： ./ieu-a-302.txt
# 2023/11/17 15:13:12 输出smr文件路径： ./ieu-a-302.ma
# 2023/11/17 15:13:18 处理时间： 5.6571224s

UKB数据

UKB数据转换为标准数据，通过UKB下载的GWAS数据与variants.tsv.bgz 通过variant列关联即可匹配出相关数据，由于消耗内存大（20GB内存以上），且时间慢的情况下。

以下为官方文件说明：

This file contains annotations on each variant in the GWAS, calculated across the analysis subset of 361,194 samples.
NOTE: The order of variants in this file matches the order of variants in the results files described below.
如何进行GWAS数据和variants.tsv.bgz的匹配？
To join these annotations with a results file, either match on the “variant” field or simply paste the columns together (e.g. “paste variants.tsv K50.gwas.imputed_v3.both_sexes.tsv”).

SNP注释文件👉：variants.tsv.bgz

示例处理 100003_irnt.gwas.imputed_v3.both_sexes.tsv.bgz 源文件大小是512.5MB


# 加载包
library(DrugTargetMR)

# 调用函数
DrugTargetMR::prepare_ukb(
            gwas_path = "./100003_irnt.gwas.imputed_v3.both_sexes.tsv.bgz",
            variants_path = "./variants.tsv.bgz",
            out_path="./",
            generate_mr = T,
            generate_smr = T)

# 2023/11/17 15:26:36 ukb gwas数据文件路径： ./100003_irnt.gwas.imputed_v3.both_sexes.tsv.bgz
# 2023/11/17 15:26:36 ukb variants数据文件路径： ./variants.tsv.bgz
# 2023/11/17 15:26:36 输出mr文件路径： ./100003_irnt.gwas.imputed_v3.both_sexes.tsv.txt
# 2023/11/17 15:26:36 输出smr文件路径： ./100003_irnt.gwas.imputed_v3.both_sexes.tsv.ma
# 2023/11/17 15:27:49 处理时间： 1m13.1463221s

FinnGen数据

示例处理finngen_R9_AB1_ACTINOMYCOSIS.gz，源文件的大小是679MB。

# 加载包
library(DrugTargetMR)

# 调用函数
DrugTargetMR::prepare_finngen(file_path = "finngen_R9_AB1_ACTINOMYCOSIS.gz",
                              out_path = "./")

# 2023/11/18 13:32:54 finngen数据文件路径： finngen_R9_AB1_ACTINOMYCOSIS.gz
# 2023/11/18 13:32:54 输出mr文件路径： finngen_R9_AB1_ACTINOMYCOSIS.txt
# 2023/11/18 13:32:54 输出smr文件路径： finngen_R9_AB1_ACTINOMYCOSIS.ma
# 2023/11/18 13:33:19 处理时间： 25.3894139s

GWAS Catalog数据

Gwas Catalog收集的数据内容，数据格式不统一，不能直接处理成MR或者SMR输出文件。转换之前，需查询数据文件列名，判断是否符合mr或smr分析要求，因此这部分需要指定需要处理的文件列名与输出文件列名的关系。固定输出列名包含：SNP, effect_allele, other_allele, eaf, beta, se, pval, samplesize, chr, pos；如果列名不指定，文件转换后的内容为NaN。

备注：来自Gwas Catalog的数据，可以使用 prepare_gwascatalog() 完成数据预处理，或者使用 prepare_others()完成数据预处理。prepare_others()适用范围更广，用法参考：其他来源数据。

示例处理28604730-GCST004744-EFO_0000571-Build37.f.tsv.gz，源文件的大小是270MB。

# step1 加载包
library(DrugTargetMR)

# step2 读取文件部分内容，判断列名映射关系

read_data <- data.table::fread("28604730-GCST004744-EFO_0000571-Build37.f.tsv.gz",
                               nrows = 10,fill = T)
colnames(read_data)
# [1] "variant_id"    "p_value"     "chromosome"  "base_pair_location"     
# [5] "odds_ratio"    "ci_lower"    "ci_upper"      "beta"                   
# [9] "standard_error" "effect_allele"  "other_allele"   "effect_allele_frequency"

# step3 调用函数
# 根据列名进行映射，不存在的列名留空; 
DrugTargetMR::prepare_gwascatalog(
  file_path = "28604730-GCST004744-EFO_0000571-Build37.f.tsv.gz",
  generate_mr = T,
  generate_smr = F,
  out_path = "./",
  col_SNP = "variant_id",
  col_effect_allele = "effect_allele",
  col_other_allele = "other_allele",
  col_eaf = "effect_allele_frequency",
  col_beta = "beta",
  col_se = "standard_error",
  col_pval = "p_value",
  col_samplesize = "",
  col_chr = "chromosome",
  col_pos = "base_pair_location"
)

# 2023/11/17 22:21:42 gwas数据文件路径： ./28604730-GCST004744-EFO_0000571-Build37.f.tsv.gz
# 2023/11/17 22:21:42 列名映射： [variant_id:SNP p_value:pval chromosome:chr base_pair_location:pos beta:beta standard_error:se effect_allele:effect_allele other_allele:other_allele effect_allele_frequency:eaf]
# 2023/11/17 22:21:42 输出mr文件路径： ./28604730-GCST004744-EFO_0000571-Build37.f.tsv.txt
# 2023/11/17 22:21:42 输出smr文件路径：./28604730-GCST004744-EFO_0000571-Build37.f.tsv.ma
# 2023/11/17 22:21:51 处理时间： 9.1701467s

其他来源数据

有些来源的数据并不能直接通过已有的预处理函数调用，因此需要根据具体的文件内容做具体的处理，转换之前，需查询数据文件列名，判断是否符合mr或smr分析要求，因此这部分需要指定需要处理的文件列名与输出文件列名的关系。固定输出列名包含：SNP, effect_allele, other_allele, eaf, beta, se, pval, samplesize, chr, pos；如果列名不指定，文件转换后的内容为NaN。

将其他来源数据快速转换为标准的mr或smr输出文件（低内存，速度快）。

注意：如果遇到输出mr或smr标准文件全是NaN，可能需要指定文件的分隔符（默认为”\t”）。

示例处理ENIGMA3_mixed_se_wo_Mean_Full_Thickness_20190429_noGC.txt.gz，源文件的大小是173MB

# step1 加载包
library(data.table)
# step2 读取文件部分内容，判断列名映射关系
read_data <- data.table::fread("ADNI_hippo_GWAS_EduCovar_NewUnrelated.assoc.linear.gz",
                               nrows = 10)
colnames(read_data)
 # [1] "SNP"    "A1"     "A2"     "FREQ1"  "BETA"   "SE"     "P"      "N"     
 # [9] "MARKER" "CHR"    "BP"    

# step3 调用函数
# 根据列名进行映射，不存在的列名留空; 
prepare_others(
  file_path = "ENIGMA3_mixed_se_wo_Mean_Full_Thickness_20190429_noGC.txt.gz",
  # 此文件分隔符为" "(默认分隔符为"\t")
  sep = " ",
  col_SNP = "SNP",
  col_effect_allele = "A1",
  col_other_allele = "A2",
  col_eaf = "FREQ1",
  col_beta = "BETA",
  col_se = "SE",
  col_pval = "P",
  col_samplesize = "N",
  col_chr = "CHR",
  col_pos = "BP"
)

# 2023/12/30 15:55:08 其他来源数据文件路径： ENIGMA3_mixed_se_wo_Mean_Full_Thickness_20190429_noGC.txt.gz
# 2023/12/30 15:55:08 列名映射： [SNP:SNP A1:effect_allele A2:other_allele FREQ1:eaf BETA:beta SE:se P:pval N:samplesize CHR:chr BP:pos]
# 2023/12/30 15:55:08 输出mr文件路径： ./ENIGMA3_mixed_se_wo_Mean_Full_Thickness_20190429_noGC.txt
# 2023/12/30 15:55:08 输出smr文件路径： ./ENIGMA3_mixed_se_wo_Mean_Full_Thickness_20190429_noGC.ma
# 2023/12/30 15:55:28 处理时间： 20.3737883s

deCODE2021 血浆蛋白

示例处理15506_34_LRP12_LRP12.txt.gz，源文件的大小是908MB。索引文件assocvariants.annotated.index，文件的大小是4GB。

# 加载包
library(DrugTargetMR)

# 调用函数
DrugTargetMR::prepare_decode(file_path = "15506_34_LRP12_LRP12.txt.gz",
                           snp_annotation_path = "assocvariants.annotated.index",
                            out_path="./",
                            generate_mr = T,
                            generate_smr = T)

# 2023/12/13 23:32:21 decode数据文件路径： 15506_34_LRP12_LRP12.txt.gz
# 2023/12/13 23:32:21 输出mr文件路径： ./15506_34_LRP12_LRP12.txt
# 2023/12/13 23:32:21 输出smr文件路径： ./15506_34_LRP12_LRP12.ma
# 2023/12/13 23:33:22 处理时间： 1m0.6360272s

assocvariants.annotated.index：由官方提供的assocvariants.annotated.txt.gz进行二次处理（包含删除错误信息，重新匹配正确数据等）；并做成索引，加快文件解析数据、降低内存消耗。

UKB-PPP2023 血浆蛋白

示例处理AAMDC_Q9H7C9_OID30236_v1_Cardiometabolic_II.tar，源文件的大小是520MB。索引文件ukp_ppp.index，文件的大小是6GB。

# 加载包
library(DrugTargetMR)

# 调用函数
# is_pos19 是否输出pos19数据，默认F
DrugTargetMR::prepare_ukb_ppp(file_path = "./AAMDC_Q9H7C9_OID30236_v1_Cardiometabolic_II.tar",
                              snp_annotation_path = "./ukp_ppp.index",
                              is_pos19 = F,
                              out_path = "./")

# 2023/12/20 23:37:00 UkbPPP数据文件路径： ./AAMDC_Q9H7C9_OID30236_v1_Cardiometabolic_II.tar
# 2023/12/20 23:37:00 输出mr文件路径： ./index/AAMDC_Q9H7C9_OID30236_v1_Cardiometabolic_II_h38.txt
# 2023/12/20 23:37:00 输出smr文件路径： ./index/AAMDC_Q9H7C9_OID30236_v1_Cardiometabolic_II_h38.ma
# 2023/12/20 23:38:12 处理时间： 1m12.5709083s

ukp_ppp.index：由官方提供的SNP RSID maps进行二次处理,并做成索引，加快文件解析数据、降低内存消耗。

标准数据处理