ALPACA on a protected research environment

[seedgeorge]

Hi,
I'm quite keen to run ALPACA, as I have already enjoyed the results of CONIPHER and RefPhase. However I am using the Genomics England RE and am not able to use Gurobi due to licensing issues (confirmed by GEL and Gurobi). Is there any scope for using an alternative solver in the future?
Kind regards,
George

[wlippa]

Hi George, and thank you for your interest in the tool and for raising this important issue. I will try to add support for alternative, non-licensed solvers and let you know here once I manage to do so.

[seedgeorge]

Thank you very much for looking into this - I appreciate it may be a bit of a hassle...

[wlippa]

Hi George, I added support for open-source solvers via Pyomo. Please try the most recent version. I suggest you use GLPK solver, as it produced slightly better results compared to SCIP, but the difference was very small. Both these solvers failed on a small number of simulated copy number segments, and I currently don't know why. In the simulations, this was a very small number (14 and 490 for GLPK and SCIP, respectively, out of a total of 90220 segments tested). Please let me know if you have any issues. The README contains the instructions describing how to use a different solver.

[seedgeorge]

Thank you for implementing this! I'm experimenting now.
My method for getting ALPACA into the research environment has been to:

• Make a docker image based on this repo
• Use docker hub to import it into GEL and run it with singularity
• Update some R packages (worth double checking that in the build instructions, I had a crash because GenomicRanges needed to be updated or something)

I've tried the input_conversion step first, and I appreciate that doesn't even include the solver issue, but unfortunately I have a crash at the 'converting refphase' step.

Converting REFPHASE output
/opt/conda/envs/alpaca/lib/python3.13/site-packages/alpaca/scripts/submodules/alpaca_input_formatting/convert_refphase_output/convert_refphase.py:91: DtypeWarning: Columns (0: seqnames) have mixed types. Specify dtype option on import or set low_memory=False.
  refphase_snps = pd.read_csv(args.refphase_snps, sep="\t")
Calculating confidence intervals for cORP00000001
Traceback (most recent call last):
  File "/opt/conda/envs/alpaca/lib/python3.13/site-packages/alpaca/scripts/submodules/alpaca_input_formatting/convert_refphase_output/convert_refphase.py", line 172, in <module>
    confidence_intervals['chr'] = confidence_intervals['segment'].apply(lambda x: int(x.split('_')[0]))
  File "/opt/conda/envs/alpaca/lib/python3.13/site-packages/pandas/core/series.py", line 5084, in apply
    ).apply()
  File "/opt/conda/envs/alpaca/lib/python3.13/site-packages/pandas/core/apply.py", line 1520, in apply
    return self.apply_standard()
  File "/opt/conda/envs/alpaca/lib/python3.13/site-packages/pandas/core/apply.py", line 1578, in apply_standard
    mapped = obj._map_values(mapper=curried)
  File "/opt/conda/envs/alpaca/lib/python3.13/site-packages/pandas/core/base.py", line 1020, in _map_values
    return arr.map(mapper, na_action=na_action)
  File "/opt/conda/envs/alpaca/lib/python3.13/site-packages/pandas/core/arrays/base.py", line 2692, in map
    return map_array(self, mapper, na_action=na_action)
  File "/opt/conda/envs/alpaca/lib/python3.13/site-packages/pandas/core/algorithms.py", line 1710, in map_array
    return lib.map_infer(values, mapper)
  File "pandas/_libs/lib.pyx", line 3071, in pandas._libs.lib.map_infer
  File "/opt/conda/envs/alpaca/lib/python3.13/site-packages/alpaca/scripts/submodules/alpaca_input_formatting/convert_refphase_output/convert_refphase.py", line 172, in <lambda>
    confidence_intervals['chr'] = confidence_intervals['segment'].apply(lambda x: int(x.split('_')[0]))
ValueError: invalid literal for int() with base 10: 'X'

I think it's because I have chr prefix for refphase, but not conipher. I'll rerun refphase, stripping chr prefix, and see how we go?
UPDATE: This doesn't change it.

I note in the example data that the refphase output doesn't have sex chromosome copy number data. I will try refphase, stripping these out.
UPDATE: This does work.

[wlippa]

Hi George, thanks for the feedback - it was indeed the X chromosome (not only the 'chr' prefix). I pushed a fix to sanitise the chromosome names, and it should now work with various conventions. 'X' will be converted to int(23) and 'Y' into int(24). Unfortunately, you will need to either pull (alpaca-0.2.6.dev7+gbef0c6522.d20260126-py3-none-any.whl) and recreate your environment or rename the chromosomes in your input data.

[wlippa]

ahh, I read the last of your updates now, happy you managed to run this part and apologies for the difficulties!

[seedgeorge]

All good, thanks for pushing a fix.

I'm testing my happily-converted data on a full apaca run now, using the parameters:
--solver pyomo
--pyomo_solver glpk
It's sitting at 'Processing Files 0%' but hasn't fallen over, and seems to have written files into the output_directory.

Some more questions:

• can I use the same directory for input and output?
• how much memory/cpu does alpaca require?
• how long should I expect the run to take?

[wlippa]

Great to hear that it worked!

1. Yes, you can use the directory for input and output.
2. Since it operates segment-by-segment, the requirements are low: I usually run it with 4GB and 2cpus, occasionally bumping to 8GB.
3. From my experience, time mostly scales with (see Extended Data Figure 2 in the paper https://www.nature.com/articles/s41586-025-09398-w/figures/8) with the number of clones and number of segments. If you want to get an estimate, find a representative tumour, keep just a few segments and run a benchmark. If your cohort is large and/or complex, or if you anticipate re-running it multiple times, you could try parallelising the runs. There are two ways to do it. A) Using only the code from this repository, run a separate instance for each tumour. Each instance will iterate over each segment. In this approach, you will need to wait for your most complex tumour to have all the results. B) If you are comfortable with Nextflow, I have a separate repo (https://github.com/McGranahanLab/ALPACA-nextflow) with an efficient, but somewhat complex, NF pipeline to parallelise over segments. The documentation there might be a bit out of date, but if you are interested, I can help you with setting this up.

[seedgeorge]

Hi, thanks for this - very helpful. Understood on all points!

For full disclosure, I'm running this from conipher and refphase-created inputs (that I ran from the outputs of the hartwig oncoanalyser pipeline - mutation calls and allele-specific copy number and so on). One quick query on that front: for refphase, should the input be rounded or unrounded major/minor copy number segments?

[wlippa]

Hi George, you should mostly use unrounded, fractional copy numbers. If you round a copy number for a segment, you are losing the information contained in the fractional part of the number and equating segments where the fractional copy number is 1.5 with segments where the fractional copy number is 1.1. The uncertainty in your fractional copy-number estimates should be captured by the confidence intervals, and you should probably check if they conform to your expectations. For example, a segment with a fractional copy number of 1.1 with many SNPs and with very narrow LogR should have smaller confidence intervals compared to another segment with 1.1 with just a handful of SNPs and high LogR spread. You can modify your CI table at will (after you get it from the pre-processing), and you can also use '--min_ci 0.5' argument to ensure that intervals are not too narrow. If you have extremely small intervals around your fractional value, say lower_ci = 1.01 and upper_ci=1.03 you will overfit to noise.

[seedgeorge]

Understood! I meant for the input to refphase, but I hear your warnings about messing with the inputs to alpaca.

For the sanitisation of the chromosome inputs, I take it if I pull the master branch now that will be resolved? I'm tweaking my docker image at the moment.

EDIT: I'm trying the current branch of ALPACA and running into a new issue... see screenshot.

[wlippa]

Hi George, re chromosome names, these should be sanitised now, i.e. we can handle string prefixes (e.g. 'chr') and X/Y/MT chromosomes (converted into 23/24/25 etc). Regarding the error, could you check if the failed listed segments are a subset of your segments, or if we are failing with all the segments? Looking at chromosome numbers, it seems like the latter might be the case. I also added detailed logs, so you could try pulling from the 'debug_reports' branch and adding '--debug' to the 'input-conversion' command and sharing the output here. I will have a look and hopefully spot the issue. Reports are written in the same directory as the output under ''convert_refphase_debug_report.txt". My hypothesis is that a merge function fails somewhere and we end up with an empty dataframe (hence why so many "segments are missing in the output").