Multivariate Fine-Mapping for multiple genes#

Description#

Multi gene fine-mapping and TWAS may also be conducted with our pipeline. This considers multiple genes jointly within specific TAD windows.

This step is similar to the multivariate fine-mapping with two main differences. 1) TAD windows with multiple genes need to be defined. The --pheno_id_map_file parameter is used for this. 2) To speed things up, the genes are filtered out if they don’t have a univariate fine mapped region. Genes may also be filtered out if they do have a univariate fine-mapped signal, but the signal is nowhere close to that of other genes. The --skip-analysis-pip-cutoff parameter is used for this.

Input#

--genoFile: path to a plink bed file containin genotypes. Include the .bed

--phenoFile: a tab delimited file containing chr, start, end, ID and path for the regions. For example:

#chr    start   end     ID      path
chr12   389319  389320  ENSG00000073614 $PATH/snuc_pseudo_bulk.Mic.mega.normalized.log2cpm.bed.gz
chr12   752578  752579  ENSG00000060237 $PATH/snuc_pseudo_bulk.Mic.mega.normalized.log2cpm.bed.gz

--covFile: path to a gzipped file containing covariates in the rows, and sample ids in the columns.

--customized-association-windows: a tab delimited file containing chr, start, end, and ID regions. For example:

#chr    start   end     TAD_id
chr1    0       10985501        chr1_0_10985501
chr1    5101787 11630758        chr1_5101787_11630758

--phenotype-names: the names of the phenotypes if multiple are included. There should be one for each phenotype file you include.

--max-cv-variants: maximum number of variants for cross-validation.

--ld_reference_meta_file: path to file containing chrom, start, end and path for linkage disequilibrium region information. For example:

#chrom  start   end     path
chr1    101384274       104443097       chr1/chr1_101384274_104443097.cor.xz,chr1/chr1_101384274_104443097.cor.xz.bim
chr1    104443097       106225286       chr1/chr1_104443097_106225286.cor.xz,chr1/chr1_104443097_106225286.cor.xz.bim

--independent_variant_list: a gzipped file containing variant information. These should be independent from one another in terms of linkage disequilibrium. For example:

chrom   pos     alt     ref     variant_id
chr1    16206   T       A       chr1:16206:T:A
chr1    16433   C       G       chr1:16433:C:G

--fine_mapping_meta: A file containg a list of gene and region information and other conditions. For example:

#chr    start   end     region_id       TSS     original_data   combined_data   combined_data_sumstats  conditions      conditions_top_loci
chr1    0       6480000 ENSG00000008128 1724356 KNIGHT_pQTL.ENSG00000008128.univariate_susie_twas_weights.rds,MiGA_eQTL.ENSG00000008128.univariate_susie_twas_weights.rds,MSBB_eQTL.ENSG00000008128.univariate_susie_twas_weights.rds,ROSMAP_Bennett_Klein_pQTL.ENSG00000008128.univariate_susie_twas_weights.rds,ROSMAP_DeJager_eQTL.ENSG00000008128.univariate_susie_twas_weights.rds,ROSMAP_Kellis_eQTL.ENSG00000008128.univariate_susie_twas_weights.rds,ROSMAP_mega_eQTL.ENSG00000008128.univariate_susie_twas_weights.rds,STARNET_eQTL.ENSG00000008128.univariate_susie_twas_weights.rds  $PATH/Fungen_xQTL.ENSG00000008128.cis_results_db.export.rds        $PATH/Fungen_xQTL.ENSG00000008128.cis_results_db.export_sumstats.rds       Knight_eQTL_brain,MiGA_GFM_eQTL,MiGA_GTS_eQTL,MiGA_SVZ_eQTL,MiGA_THA_eQTL,BM_10_MSBB_eQTL,BM_22_MSBB_eQTL,BM_36_MSBB_eQTL,BM_44_MSBB_eQTL,monocyte_ROSMAP_eQTL,Mic_DeJager_eQTL,Ast_DeJager_eQTL,Oli_DeJager_eQTL,Exc_DeJager_eQTL,Inh_DeJager_eQTL,DLPFC_DeJager_eQTL,PCC_DeJager_eQTL,AC_DeJager_eQTL,Mic_Kellis_eQTL,Ast_Kellis_eQTL,Oli_Kellis_eQTL,OPC_Kellis_eQTL,Exc_Kellis_eQTL,Inh_Kellis_eQTL,Ast_mega_eQTL,Exc_mega_eQTL,Inh_mega_eQTL,Oli_mega_eQTL,STARNET_eQTL_Mac       Knight_eQTL_brain,MiGA_GFM_eQTL,MiGA_GTS_eQTL,MiGA_SVZ_eQTL,MiGA_THA_eQTL,BM_10_MSBB_eQTL,BM_22_MSBB_eQTL,BM_36_MSBB_eQTL,BM_44_MSBB_eQTL,monocyte_ROSMAP_eQTL,Mic_DeJager_eQTL,Ast_DeJager_eQTL,Oli_DeJager_eQTL,Exc_DeJager_eQTL,Inh_DeJager_eQTL,DLPFC_DeJager_eQTL,PCC_DeJager_eQTL,AC_DeJager_eQTL,Mic_Kellis_eQTL,Ast_Kellis_eQTL,Oli_Kellis_eQTL,OPC_Kellis_eQTL,Exc_Kellis_eQTL,Inh_Kellis_eQTL,Ast_mega_eQTL,Exc_mega_eQTL,Inh_mega_eQTL,Oli_mega_eQTL,STARNET_eQTL_Mac

--phenoIDFile: A bed file containing a list of genes and their LD region. For example:

TAD_id  ID
chr19_0_13957223        ENSG00000172270
chr19_0_13957223        ENSG00000099864
chr19_0_13957223        ENSG00000011304

--skip-analysis-pip-cutoff: A number of the pip cutoff.

--coverage

--maf

--pheno_id_map_file: A file containing IDs and genes. For example:

ID      gene
chr20:50940933:50941105:clu_44490_-:ENSG00000000419     ENSG00000000419
chr20:50940933:50941129:clu_44490_-:ENSG00000000419     ENSG00000000419
chr20:50936262:50942031:clu_44490_-:ENSG00000000419     ENSG00000000419

--prior-canonical-matrices

--save-data: whether to save intermediate data or not

--twas-cv-folds: Perform K folds valiation CV for TWAS. Set this to zero to skip

--trans-analysis: Include this if doing trans-analysis (not using phenotypic coordinate information)

--region-name: if you only wish to analyze one region, then include the ID of a region found in the customized-association-windows file

--cwd: output file path

Output#

Minimal Working Example Steps#

ii. Run the Fine-Mapping with mvSuSiE#

sos run pipeline/mnm_regression.ipynb mnm_genes \
    --name ROSMAP_mega_eQTL --cwd output/mnm_regression/mnm_genes \
    --genoFile output/genotype_by_chrom/wgs.merged.plink_qc.genotype_by_chrom_files.txt \
    --phenoFile output/phenotype/phenotype_by_chrom_for_cis/bulk_rnaseq.phenotype_by_chrom_files.region_list.txt \
    --covFile output/covariate/bulk_rnaseq_tpm_matrix.low_expression_filtered.outlier_removed.tmm.expression.covariates.wgs.merged.plink_qc.plink_qc.prune.pca.Marchenko_PC.gz \
    --customized-association-windows reference_data/TAD/TADB_enhanced_cis.bed \
    --region-name ENSG00000049246 ENSG00000054116 ENSG00000116678 ENSG00000073921 \
    --phenotype-names test_pheno \
    --max_cv_variants 5000 --skip-analysis-pip-cutoff 0.025  \
	--ld_reference_meta_file data/ld_meta_file_with_bim.tsv \
    --pheno-id-map-file data/mnm_regression/pheno_id_map_file.tsv \
    --fine-mapping-meta data/mnm_regression/fine_mapping_meta.tsv

Anticipated Results#

For each gene and region, multivariate multigene finemapping will produce a file containing results for the top hits and a file containing twas weights produced by susie.

ROSMAP_Ast_DeJager_eQTL.chr11_77324757_86627922.multigene_bvrs.rds:

  • for each region name, includes:

  1. mrmash_fitted

  2. reweighted_mixture_prior

  3. reweighted_mixture_prior_cv

  4. mvsusie_fitted

  5. variant_names

  6. analysis_script

  7. other_quantitites

  8. analysis_script

  9. top_loci

  10. susie_result_trimmed

  11. total_time_elapsed

  12. region_info

ROSMAP_Ast_DeJager_eQTL.chr11_77324757_86627922.multigene_data.rds:(from the –save-data argument)

ROSMAP_Ast_DeJager_eQTL.chr11_77324757_86627922.multigene_twas_weights.rds:

  • for each region name and for each gene within that region, includes:

  1. twas_weights - from mrmash and mvsusie

  2. twas_predictions - from mrmash and mvsusie

  3. variant_names

  4. total_time_elapsed

  5. region_info