Skip to content

Export fSuSiE curves and bands in the canonical fine-mapping workflow - #1426

Open
hsun3163 wants to merge 2 commits into
StatFunGen:mainfrom
hsun3163:fix/fsusie-functional-export
Open

hsun3163 wants to merge 2 commits into
StatFunGen:mainfrom
hsun3163:fix/fsusie-functional-export

Conversation

@hsun3163

@hsun3163 hsun3163 commented Sep 23, 2026 •

Copy link
Copy Markdown
Contributor

Regional fSuSiE results need both reusable fine-mapping summaries for colocalization and the fitted functional effects for plotting. The existing workflow only declares the regional RDS, and trimmed fSuSiE fits have discarded the curve fields needed for an effect export. This PR adds the combined top-loci/effect table to the existing workflow, paired with pecotmr #600, which retains those fields and stores each joint fit once.

fine_mapping.R

After writing a regional fSuSiE result, export the retained grid curves, native uncertainty bands, probe names/positions, and interpolated probe effects alongside credible-set variants. The table follows the existing 19-column bulk epigenetic export schema and appends grid_band_halfwidth, a semicolon-packed vector aligned with grid_positions. The saved native bands can be reconstructed as grid_effects -/+ grid_band_halfwidth; the column is not a standard error. Symmetry is checked before reducing the two band rows to half-widths.

Probe effects use linear interpolation with constant endpoint extrapolation, matching the existing epigenetic exporter. A variant belonging to multiple credible sets has one row per effect. A region with no credible sets writes the same header with no data rows. The RDS still retains all fitted effects.

The QTL region argument is converted from its CLI string into GenomicRanges::GRanges at the wrapper boundary. This uses pecotmr's supported region type and fixes the overlap-selection failure exposed by regional calls. Gene-ID and GWAS input modes do not use this branch.

The exporter selects a one-row collection and calls the public getSusieFit(), getVariantIds(), and getTopLoci() accessors. The original draft directly read x$entry[[i]], which works on 0.6.12 but returns NULL on the current GRangesList-based result. CI exposed that reader defect as getSusieFit(NULL). The fit had been saved; it was the lookup that failed.

The exporter runs only for regional fSuSiE calls and fails clearly if an old fit lacks the required functional fields. Existing trimmed files cannot regenerate discarded curves; they need refitting. Successful RDS serialization alone did not demonstrate that plotting information was available: this is an output-contract gap, not a file-format artifact.

mnm_regression.ipynb and the shared Snakemake target

fsusie_2 declares both the RDS and regional BED. fsusie_3 combines regional BEDs into one sorted, bgzipped, tabix-indexed table per configured name. It uses the same notebook and scripts as the fit; there is no separate exporter workflow.

The existing local shared fSuSiE rule module and its small phenotype-manifest builder are included here so the upstream shared Snakefile can declare the complete target, including the final table/index. This accounts for most of the added workflow lines. The rule builds the QtlDataset once per context, selects TADs, and invokes the canonical notebook. It forwards chromosome scope, memory, walltime, PC count, and post-processing. TI remains the default; the application configs now state TI explicitly. This path does not invoke TWAS-weight training.

Validation and compatibility

  • The canonical qtl_dataset_construct+fsusie SoS test passed on the committed chr22 fixture with TI and one PC per context, using both patched pecotmr 0.6.12 and patched 0.8.2. It produced four result rows (one joint fSuSiE fit and one PC fit per context), a nonempty 20-column table, aligned effects/positive band widths, and a tabix index that answered a regional query.
  • Dry-runs passed for the clean shared-Snakefile fixture and both CARDIA and MESA application configs. Their rendered commands use TI, one PC, and autosomes only.
  • The 0.6.12 check used the installed, tested pecotmr backport 7dde6d3. The exporter uses public collection accessors for both result representations.
  • Existing univariate and multivariate notebook steps retain their outputs. fSuSiE gains declared export products and therefore requires the companion pecotmr retention change. Production reruns were not launched as part of this PR.

Protocol CI installs the immutable companion pecotmr revision from PR #600 before running its existing tests. The released 0.8.2 package lacks the retained curves; testing against that release cannot validate the new export. The source pin can be replaced with a released dependency after the package change ships.

The five pipeline-level assertions in the companion package PR retain the old per-probe row expectations at the author's request. They deliberately expose the output-contract difference; this exporter repair does not change those test expectations.

This branch has not been deployed

No deployments
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

1 participant