随着处理更多的gwas数据，慢慢发现MungeSumstats包的妙处，这期就介绍一下这个包的详细参数，方便大家处理自己的数据——
代码语言：javascript复制if (!require("BiocManager", quietly = TRUE))
install.packages("BiocManager")
BiocManager::install("MungeSumstats")
第一步肯定是安装啦~
这个包配合以下几个包食用更佳：
代码语言：javascript复制BiocManager::install("SNPlocs.Hsapiens.dbSNP155.GRCh38")
BiocManager::install("BSgenome.Hsapiens.NCBI.GRCh38")
##or
BiocManager::install("SNPlocs.Hsapiens.dbSNP155.GRCh37")
BiocManager::install("BSgenome.Hsapiens.1000genomes.hs37d5")
请注意，与参考基因组一起使用的默认 dbSNP 版本是 Bioconductor 上可用的最新版本（当前为 dbSNP 155），但旧版本也可用。使用
dbSNP输入参数来控制它。
MungeSumstats 推断效应等位基因将始终是 A2 等位基因，这是IEU GWAS VCF所做的方法，并且此处也采用了这种方法。该推论首先来自输入文件的列标题，但是，等位基因翻转检查通过将 A1（应该是参考等位基因）与参考基因组进行比较来确保这一点。 如果 SNP 的 A1 DNA 碱基与参考基因组不匹配，但 A2（应该是替代等位基因）与参考基因组匹配，则等位基因将与效应信息（例如 Beta、优势比、签名汇总统计、FRQ、Z）一起翻转。 MungeSumstats可以处理 VCF、txt、tsv、csv 文件类型或这些文件类型的 .gz/.bgz 版本。 该软件包还使用户能够灵活地将重新格式化的文件导出为制表符分隔的 VCF 或 R 本机对象，例如 data.table、GRanges 或 VRanges 对象。输出还可以以LDSC 就绪格式输出，这意味着文件可以直接输入 LDSC，无需额外修改。
为什么要反复强调A1&A2，因为很多时候Allele1—— the Effect Allele；Allele2——the other Allele，这与此包相反。如果不了解这个包，可能会让我们的分析南辕北辙，回头却不知错在何处~~>_<~~
参数介绍
MungeSumstats的核心函数是format_sumstats
- convert_small_p 要将
p-values < 5e-324转换为0吗？小 p 值超过 R 限制，可能会导致 LDSC/MAGMA 出现错误，应进行转换。默认值为 TRUE。 - convert_large_p p 值 >1 是否转换为 1？P 值 >1 应该是不可能的，并且可能会导致 LDSC/MAGMA 错误，应进行转换。默认值为 TRUE。
- convert_neg_p p 值 <0 是否应该转换为 0？负 p 值不应该是可能的，并且可能会导致 LDSC/MAGMA 错误，应进行转换。默认值为 TRUE。
- compute_z 是否从 P 计算 Z 分数列。默认值为 FALSE。请注意，为每个 SNP 计算 Z 分数并不完全正确，并且可能会导致功效损失。这只能作为最后的手段。
- force_new_z 当“Z”列已经存在时，默认使用它。要从 P 设置为 TRUE 覆盖并计算新的 Z 分数列。
- compute_n 是否插补 N。默认值 0 不会插补，任何其他整数将被插补为数据集中每个 SNP 的 N（样本大小）。请注意，为每个 SNP 估算样本量并不正确，只能作为最后的手段。 N 还可以通过为该字段传递其中之一或多个向量来输入“ldsc”、“sum”、“giant”或“metal”。Sum 和整数值在输出中创建 N 列，而 Giant、metal 或 ldsc 创建 Neff 或有效样本大小。如果传递多个，则会指示用于推导它的公式。
- convert_n_int 如果N（样本数）不是整数，是否应该四舍五入？默认值为 TRUE。Analysis_trait 如果研究多个性状，则用于 GWAS 分析的性状名称。默认值为 NULL。
- impute_beta 如果sumstats中不存在BETA，是否应使用其他效果数据来估算BETA。请注意，此估算是近似值，因此可能会对下游分析产生影响。谨慎使用。MungeSumstats 将尝试估算 beta 的不同方法（按此顺序或优先级）是：1. log(OR) 2. Z x SE。默认值为 FALSE。
- es_is_beta 是否将 ES 映射到 BETA。我们将 BETA 视为任何类似 BETA 的值（包括效应大小）。如果您的 sumstats 不是这种情况，请将其更改为 FALSE。默认值为 TRUE。
- impute_se 如果 sumstats 中不存在标准误差，是否应使用其他效应数据来估算标准误差。请注意，此估算是近似值，因此可能会对下游分析产生影响。谨慎使用。MungeSumstats 将尝试估算的不同方法（按此顺序或优先级）是：1.BETA / Z；2.绝对值绝对值（BETA/qnorm（P/2））。默认值为 FALSE。
- analysis_trait 如果研究多个性状，则用于 GWAS 分析的性状名称。默认值为 NULL。
- INFO_filter 插补信息分数允许的最小值（如果在 sumstatsfile 中存在）。默认 0.9
- FRQ_filter 0-1 SNP 频率（FRQ）允许的最小值（即等位基因频率（AF））（如果在 sumstats 文件中存在）。默认情况下不进行过滤，即值为 0。
- pos_se 是否应该检查标准错误 (SE) 列以确保它大于 0？那些存在的内容将被删除（如果 sumstats 文件中存在）。默认为TRUE。
- effect_columns_nonzero 应检查数据BETA、OR（比值比）、LOG_ODDS、SIGNED_SUMSTAT 中的效果列，以确保没有 SNP=0。那些这样做的被删除（如果存在于 sumstats 文件中）。默认为TRUE。
- N_std 需要删除高于 SNP N 平均值的标准差数。默认值为 5。N_dropNA控制是否删除缺少 N 值的 SNP（默认值为 TRUE）。N_dropNA 删除缺少 N 的行。默认值为 TRUE。
- rmv_chr向量或字符 应删除 SNP 的染色体。如果不需要过滤，则使用 NULL。默认为 X、Y 和线粒体。
- rmv_chrPrefix 控制是否从染色体名称中删除“chr”/“CHR”（默认为 TRUE）。
- on_ref_genome 应检查所有 SNP 是否均按 SNP ID 位于参考基因组上。任何不在参考基因组上的 SNP 将使用染色体和碱基对位置数据从参考基因组（如果可能）进行校正。默认为 TRUE
- Convert_ref_genome要转换的参考基因组的名称（“GRCh37”或“GRCh38”）。仅当当前基因组构建不匹配时才会发生这种情况。默认不转换基因组构建（NULL）。
- strand_ambig_filter 应删除具有链模糊等位基因的 SNP。默认为FALSE。
- allele_flip_check 是否应根据参考基因组检查等位基因列以推断是否需要翻转。默认值为 TRUE。
- allele_flip_drop控制是否删除 A1 或 A2 碱基对值均不与参考基因组匹配的 SNP。默认值为 TRUE。
- allele_flip_z 控制 Z 分数值是否应与effect 和 FRQ 列（例如 Beta）一起翻转。默认值为 TRUE。
- allele_flip_frq控制频率 (FRQ) 值是否应与effect 和 Z 分数列（例如 Beta）一起翻转。默认值为 TRUE。
- bi_allelic_filter 应删除非双等位基因 SNP。默认为 TRUE
- snp_ids_are_rs_ids 如果输入的 SNP ID 被推断为 RS ID 或某个任意 ID。默认值为 TRUE。有时，汇总统计信息可以在一行上有多个 RSID（即与一个 SNP 相关），例如“rs5772025_rs397784053”。这可能会导致错误，因此默认情况下，将保留第一个 RS ID，并删除其余的，例如“rs5772025”。如果您只想完全删除这些 SNP，请将其设置为 TRUE。默认值为 FALSE。
- frq_is_maf 传统上 FRQ 列旨在显示次要/影响等位基因频率 (MAF)，但有时可以将主要等位基因频率推断为 FRQ 列。该逻辑变量指示如果频率值似乎与主要等位基因相关，即 >0.5，则 FRQ 列应重命名为 MAJOR_ALLELE_FRQ。默认情况下不会发生映射，即为 TRUE。
- indels 您的 Sumstats 文件是否包含 Indel？这些不存在于我们的参考文件中，因此如果该值为 TRUE，它们将被排除在检查之外。默认值为 TRUE。
- dbSNP 用作参考的 dbSNP 版本 - 默认为可用的最新版本 (155)。请注意，dbSNP 155 中的 SNP 比 144 中的 SNP 多了 9 倍，运行时间将会增加。
- sort_coordinates是否按结果 sumstats 的坐标排序。
- nThread用于并行进程的线程数。
- write_vcf是否写入 VCF (TRUE) 或表格文件 (FALSE)。而tabix_index是一个 输入，用于确定是否用tabix对格式化的汇总统计数据建立索引，以便快速查询。
- return_data返回
data.table，GRanges或VRanges直接返回给用户。否则，返回保存数据的路径。默认值为 FALSE。 - return_format如果 return_data 为 TRUE。要返回的对象类型（“data.table”、“vranges”、“granges”）。
- save_format 确保输出格式符合所有要求，可直接传入 LDSC（"ldsc"），而无需额外润色，或在保存为 VCF 之前传入 IEU OpenGWAS 格式（"opengwas"）。
- log_folder_ind应存储包含所有过滤掉的 SNP 的日志文件（每个过滤器单独的文件）。数据以与生成的 sumstats 文件指定的相同格式输出。
- log_mungesumstats_msgs 应该存储包含 MungeSumstats 在运行中打印的所有消息和错误的日志。
- imputation_ind 应该为每个插补步骤添加一列，以显示哪些 SNP 对不同字段具有插补值。这包括表示 SNP 等位基因翻转（翻转）的字段。对于翻转值，这表示等位基因是否根据 MungeSumstats 从输入列标题中选择的 A1、A2 进行切换，因此可能与创建者的意图不符。请注意，这些列将出现在返回的格式化摘要统计信息中。
- log_folder日志文件和要存储的 MungeSumstats 消息日志的目录路径。默认是临时目录。如果存在同名的格式化文件，则将跳过格式化并导入该文件（默认）。设置为覆盖此设置。
- mapping_file MungeSumstats 有一个预定义的列名映射文件，该文件应涵盖最常见的列标题及其解释。但是，如果 youf 文件中的列标题丢失，我们提供的映射不正确，您可以提供自己的映射文件。必须是 2 列数据框，列名称为“未更正”和“已更正”。请参阅
data(sumstatsColHeaders)默认映射和必要的格式。
基因组转换
代码语言：javascript复制MungeSumstats 将 liftover() 函数作为通用工具提供给用户。
sumstats_dt_hg38 <- MungeSumstats::liftover(sumstats_dt = sumstats_dt,
ref_genome = "hg19",
convert_ref_genome = "hg38")
文件类型转换和列标题标准化
- 读入数据并规范标题名称
dat <- MungeSumstats::read_sumstats(path = eduAttainOkbayPth,
standardise_headers = TRUE)
- 以压缩、制表符分隔、制表符索引文件的形式写入磁盘
formatted_path <- MungeSumstats::write_sumstats(sumstats_dt = dat,
save_path = formatted_path,
tabix_index = TRUE,``
write_vcf = FALSE,
return_path = TRUE)
从 Open GWAS 导入 GWAS 统计摘要
- 根据疾病搜索
metagwas <- MungeSumstats::find_sumstats(traits = c("parkinson","alzheimer"),
min_sample_size = 1000)
head(metagwas,3)
ids <- (dplyr::arrange(metagwas, nsnp))$id
- 根据ieu id
datasets <- MungeSumstats::import_sumstats(ids = "ieu-a-298",
ref_genome = "GRCH37")
一点小建议：虽然官方指南说这个包可以很好的处理vcf文件，但我尝试了一下，不建议这么做，很慢且极耗内存……
参考：
https://github.com/neurogenomics/MungeSumstats
https://neurogenomics.github.io/MungeSumstats/articles/MungeSumstats.html
