TWAS and COLOC integration (INTACT) - FunGen-xQTL Consortium

Run INTACT on the PTWAS and fastenloc output to integrate transcriptome-wide association (TWAS) evidence with colocalization (COLOC) evidence into a single gene-level posterior probability.

Overview¶

The posterior probabilities in the INTACT function are calculated through the following steps:

Compute Bayes Factors: The function converts TWAS z-scores into Bayes factors by computing a grid of Bayes factors using different values of K (a vector of values over which Bayesian model averaging is performed, default c(1,2,4,8,16)), then performing Bayesian model averaging using the log sum exp trick.
Compute Prior Probabilities: The prior_fun function converts GLCP into prior probabilities. prior_fun is set as linear by default.
Compute Posterior Probabilities.

When to run it. After PTWAS and fastenloc have both been run on the same gene set, to combine their evidence.

Input¶

--ptwas-file: PTWAS output, a TSV file with (at minimum) columns GENE, STAT (TWAS z-score) and SUBCLASS.
--fastenloc-file: fastenloc output, a whitespace-delimited table with columns Gene and GLCP (gene-level colocalization probability).
--tissue: tissue / dataset label used to name the output (e.g. DLPFC).
--alpha: FDR significance threshold (default 0.05).

The genes are matched between the two inputs by joining on the gene name.

--cwd: the directory outputs are written to.

Output¶

{cwd}/{tissue}.INTACT.rds - the merged PTWAS and fastenloc table with the INTACT posterior probability added, sorted by descending posterior. Example output/intact/DLPFC.INTACT.rds, a data frame of 40 genes and 5 columns:
```
'data.frame':  40 obs. of  5 variables:
 $ gene      : chr  "GENE021" "GENE023" "GENE004" "GENE034" ...
 $ zscore    : num  -6.55 -4.45 4.03 -3.41 -3.61 ...
 $ GLCP      : num  0.9975 0.0909 0.2677 0.8608 0.3994 ...
 $ intact_pip: num  1 0.983 0.975 0.963 0.941 ...
 $ fdr_sig   : logi  TRUE TRUE TRUE TRUE TRUE FALSE ...
```
zscore comes from the PTWAS input and GLCP from fastenloc; intact_pip is the posterior this module computes from the two, and fdr_sig flags the genes passing --alpha.

Minimal Working Example¶

One workflow, intact, run once per tissue. It requires the R package INTACT from Bioconductor, provided by the pipeline environment as bioconductor-intact.

Timing: ~1-2 min (on toy dataset)

sos run pipeline/intact.ipynb intact \
    --fastenloc-file tests/fixtures/intact/protocol_example.fastenloc.gene.out \
    --ptwas-file tests/fixtures/intact/protocol_example.ptwas.output \
    --tissue DLPFC \
    --cwd output/intact

Environment note. This step requires the R package INTACT (Bioconductor), provided by the pipeline environment as bioconductor-intact. Supply real PTWAS and fastenloc output files for a production run.

Command Interface¶

sos run pipeline/intact.ipynb -h

usage: sos run code/SoS/pecotmr_integration/intact.ipynb
               [workflow_name | -t targets] [options] [workflow_options]
  workflow_name:        Single or combined workflows defined in this script
  targets:              One or more targets to generate
  options:              Single-hyphen sos parameters (see "sos run -h" for details)
  workflow_options:     Double-hyphen workflow-specific parameters

Workflows:
  intact

Global Workflow Options:
  --cwd output (as path)
                        Workdir
  --modular-script-dir code/script (as path)
                        Directory of the modular analysis scripts (code/script)
  --fastenloc-file ''
                        fastenloc output
  --ptwas-file ''
                        ptwas output
  --tissue ''
                        dataset
  --alpha 0.05 (as float)
  --QTL eQTL
                        QTL data type
  --container ''
  --job-size 1 (as int)
  --walltime 5h
  --mem 8G
  --numThreads 1 (as int)

Sections
  intact:

Workflow implementation¶

[global]
# Workdir
parameter: cwd = path("output")
# Directory of the modular analysis scripts (code/script)
parameter: modular_script_dir = path('code/script')
# fastenloc output
parameter: fastenloc_file = ""
# ptwas output
parameter: ptwas_file = ""
# dataset 
parameter: tissue = ''
parameter: alpha = 0.05
# QTL data type
parameter: QTL = 'eQTL'
parameter: container = ''
parameter: job_size = 1
parameter: walltime = "5h"
parameter: mem = "8G"
parameter: numThreads = 1

[intact]
input: ptwas_file, fastenloc_file
output: f'{cwd}/{tissue}.INTACT.rds'
task: trunk_workers = 1, trunk_size = job_size, walltime = walltime, mem = mem, cores = numThreads, tags = f'{step_name}_{_output:bn}'
bash: expand = "${ }", stderr = f'{_output:nn}.stderr', stdout = f'{_output:nn}.stdout', container = container
    Rscript ${modular_script_dir}/pecotmr_integration/intact.R \
        --ptwas-file ${_input[0]} \
        --fastenloc-file ${_input[1]} \
        --alpha ${alpha} \
        --output ${_output}