Step 25: Polygenic Risk Scores (PRS)¶
What This Does¶
Scores you for 15 common conditions and traits with the PGS Catalog's own calculator, pgsc_calc. Each polygenic score adds up the small effects of hundreds to millions of variants into one number. With the ancestry reference panel installed, each score is also given as a percentile: where your score falls among the reference samples whose genetic ancestry is most similar to yours. Without the panel the step reports the raw sum and says "raw score only", because a raw sum from one person cannot be compared with anyone. With the panel, a score that a published study evaluated in your reference group also shows that study's effect: the odds or hazard ratio it measured, with the PGS Catalog id of the evaluation.
Why¶
Most common diseases (heart disease, diabetes, cancer) are not caused by a single gene. They come from the combined effect of many variants, each adding a small amount of risk. A score sums those contributions with weights from large genome-wide association studies (GWAS). The sum only means something next to the sums of other people scored the same way, and the distribution of sums differs between ancestries. That is why the percentile is taken among the most similar reference group, and why pgsc_calc is used: it matches variants the way the PGS Catalog intends, and with a reference panel it reports ancestry-adjusted percentiles.
Tool¶
- pgsc_calc (Lambert et al., Nat Genet 2024), the PGS Catalog Calculator, release
PGSC_CALC_VERSIONinversions.env, Apache-2.0. It is a Nextflow pipeline of its own; step 25 starts it. - PGS Catalog: curated repository of published polygenic scoring files.
Docker Image¶
pgsc_calc runs its steps in these images, all pinned in versions.env: PGSC_UTILS_IMAGE (matching and ancestry adjustment), PLINK2_IMAGE (scoring), PGSC_FRAPOSA_IMAGE (projection onto the panel), PGSC_ZSTD_IMAGE, PGSC_PYYAML_IMAGE and PGSC_REPORT_IMAGE (its HTML report). Step 25 also uses PYTHON_IMAGE and BCFTOOLS_IMAGE. Image versions lists the tags.
Input¶
- VCF from DeepVariant (step 3):
${GENOME_DIR}/${SAMPLE}/vcf/${SAMPLE}.vcf.gz - gVCF from DeepVariant (step 3), when present:
${GENOME_DIR}/${SAMPLE}/vcf/${SAMPLE}.g.vcf.gz. The score positions are genotyped from it, so sites where you match the reference count. - The scores of
assets/pgs_scores.tsv, downloaded bysetup.shinto${GENOME_DIR}/prs_scores/. - Optional: the ancestry reference panel,
${GENOME_DIR}/reference/pgsc_calc/pgsc_1000G_v1.tar.zstwith its site list beside it (scripts/setup.sh --ancestry-panel <genome_dir>): 7.4 GB, and about 23 GB of disk while pgsc_calc uses it (measured below). - Java 17+ and Nextflow on the host, as for
run-all.sh.
Command¶
./scripts/setup.sh /path/to/genome_dir # the scores, pgsc_calc and its plugin (once)
./scripts/setup.sh --ancestry-panel /path/to/genome_dir # optional: percentiles (about 7 GB)
./scripts/25-prs.sh your_name male # or female: the scores' chrX rows are scored too
./scripts/25-prs.sh your_name # no sex: chrX rows are left out
PRS_KEEP_AMBIGUOUS=true ./scripts/25-prs.sh your_name male # opt-in: keep the A/T and C/G rows
ANCESTRY_PANEL=none scores without the panel even when it is installed. PGSC_MAX_MEMORY (for example 24.GB) caps the memory each pgsc_calc task may use; the default is three quarters of the machine's RAM. The default scores need about 14 to 18 GB of it, see Memory.
PRS_KEEP_AMBIGUOUS=true (Nextflow: --prs_keep_ambiguous) keeps each score's strand-ambiguous rows, see Strand-ambiguous rows. The default is false. It needs step 3's gVCF: the step refuses it on a sample with only the variant VCF, and so does the PRS process.
Conditions Scored¶
The list is assets/pgs_scores.tsv. Each label is the trait_reported value the PGS Catalog REST API returns for that score; scripts/ci/check-pgs-labels.sh compares the two, so an ID cannot be printed under the wrong disease.
| Condition (PGS Catalog trait) | PGS ID | Variants | Publication |
|---|---|---|---|
| Coronary artery disease | PGS000018 | 1,745,179 | Inouye et al. 2018, J Am Coll Cardiol |
| Type 2 diabetes (T2D) | PGS000014 | 6,917,436 | Khera et al. 2018, Nat Genet |
| Breast cancer | PGS000004 | 313 | Mavaddat et al. 2018, Am J Hum Genet |
| Prostate cancer | PGS000662 | 269 | Conti et al. 2021, Nat Genet |
| Atrial fibrillation | PGS000016 | 6,730,541 | Khera et al. 2018, Nat Genet |
| Late-onset Alzheimer’s disease | PGS000334 | 22 | Zhang et al. 2020, Nat Commun |
| Body mass index (BMI) | PGS000027 | 2,100,302 | Khera et al. 2019, Cell |
| Inflammatory bowel disease | PGS000017 | 6,907,112 | Khera et al. 2018, Nat Genet |
| Colorectal cancer | PGS000055 | 76 | Schmit et al. 2019, J Natl Cancer Inst |
| Systolic blood pressure | PGS002257 | 884 | Evangelou et al. 2018, Nat Genet |
| LDL cholesterol | PGS000065 | 103 | Kuchenbaecker et al. 2019, Nat Commun |
| HDL cholesterol | PGS000064 | 120 | Kuchenbaecker et al. 2019, Nat Commun |
| Triglycerides (TG) | PGS000066 | 101 | Kuchenbaecker et al. 2019, Nat Commun |
| Asthma | PGS002727 | 985,837 | Namjou et al. 2022, J Allergy Clin Immunol |
| Height | PGS002802 | 1,103,042 | Yengo et al. 2022, Nature |
The blood pressure score comes from a GWAS of over a million people and has been evaluated by four other publications in the catalog. The three lipid scores come from a study of how lipid loci transfer across ancestries: they were evaluated in five cohorts of three ancestries. All four are GRCh38-harmonised, additive and on the autosomes only.
The asthma score is the one the eMERGE network chose for its asthma result (Lennon et al. 2024, Nat Med, PMID 38374346, Supplementary Table 3). Height is there as a check of the engine on a trait everyone has and that is easy to measure; its evaluations in the catalog report no odds or hazard ratio, so it shows no effect. eMERGE's chronic kidney disease and type 1 diabetes scores are not in the list. Its kidney model is the Khan et al. 2022 score (PGS002237) plus a recessive APOL1 term, so it cannot be scored as published from additive weights alone. Its type 1 diabetes score starts from PGS000024, whose file holds HLA interaction terms. This pipeline scores additive weights only, so it cannot score either one.
The third column of the file, applies_to, names the sex a score was built for: female for PGS000004 (breast cancer), male for PGS000662 (prostate cancer), empty for the rest. When you declare the other sex, the score keeps its raw sum and gets no percentile, and its Sex_Check says not applicable to the declared sex. A percentile of breast cancer for a man, or prostate cancer for a woman, would read as a risk the score was never built to give. Without a declared sex the two scores are reported as usual, with Sex_Check no sex declared and a line in the reports.
Effect sizes and the ancestry they were measured in¶
A polygenic score is built and tested mostly in one ancestry, and its effect in other ancestries can be smaller or larger. The PGS Catalog records each published evaluation of a score as a performance metric with its own id (PPM), the samples it was measured in and the effect it found. The last three columns of the file name the evaluations the reports use, one per ancestry, in the same order in each column:
ppm_id: the evaluations, for examplePPM000591;PPM000594;PPM000597for PGS000018.ancestry: the 1000 Genomes group each was measured in, as pgsc_calc names your reference group:EUR;AFR;AMR. The catalog's broad ancestry terms map to them as "European" to EUR; "African unspecified", "African American or Afro-Caribbean" and "Sub-Saharan African" to AFR; "Hispanic or Latin American" to AMR; "East Asian" to EAS; "South Asian" to SAS. An evaluation in mixed, unreported or other ancestries is not used.effect: the odds ratio (OR) or hazard ratio (HR) the catalog records for it:HR 1.53;HR 1.27;HR 1.53. It is in the unit the evaluation used. Open the PPM on the catalog (https://www.pgscatalog.org/rest/performance/PPM000591, or the score's page) to read it.
The evaluations used today:
| Score | Evaluations (ancestry, effect) | Publication |
|---|---|---|
| PGS000018, coronary artery disease | PPM000591 (EUR, HR 1.53), PPM000594 (AFR, HR 1.27), PPM000597 (AMR, HR 1.53) | Dikilitas et al. 2020, Am J Hum Genet |
| PGS000014, type 2 diabetes | PPM021437 (EUR, OR 1.71), PPM021589 (AFR, OR 1.27), PPM021615 (EAS, OR 1.28), PPM021635 (SAS, OR 1.48) | Ritchie et al. 2024, medRxiv (preprint) |
| PGS000004, breast cancer | PPM000004 (EUR, OR 1.61) | Mavaddat et al. 2018, Am J Hum Genet |
| PGS000662, prostate cancer | PPM022982 (EUR, HR 1.72), PPM015482 (AFR, OR 1.86) | Liu et al. 2024; Kim et al. 2022 |
| PGS000016, atrial fibrillation | PPM014988 (EUR, HR 1.4) | Marston et al. 2022 |
| PGS000334, late-onset Alzheimer's disease | PPM022981 (EUR, HR 1.93) | Liu et al. 2024 |
| PGS002727, asthma | PPM014751 (EUR, OR 1.57), PPM014752 (AFR, OR 1.43), PPM014753 (AMR, OR 1.63), PPM014754 (EAS, OR 1.32), PPM014755 (SAS, OR 1.32) | Namjou et al. 2022, J Allergy Clin Immunol |
The other eight have none (checked on 2026-11-03). BMI, blood pressure, the lipids and height are measurements: their scores are evaluated with a beta or a variance explained, and an odds ratio for a related diagnosis (obesity, hypertension) is not the effect on the trait the score measures, so none is used. The catalog has no OR or HR evaluation of PGS000017 (inflammatory bowel disease) or PGS000055 (colorectal cancer).
scripts/ci/check-pgs-labels.sh asks https://www.pgscatalog.org/rest/performance/
There is no schizophrenia score yet. An earlier version listed PGS000738 as schizophrenia, but that score is for vitiligo; a schizophrenia row comes back only once a score is chosen from the catalog and checked against the API.
To add a score, add a line to assets/pgs_scores.tsv with the label the API gives, female or male in the third column when the score is for one sex only, and its evaluations when it has some:
curl -s https://www.pgscatalog.org/rest/score/PGS000017 | jq -r '.trait_reported, .variants_number'
url='https://www.pgscatalog.org/rest/performance/search?pgs_id=PGS000017'
while [ -n "$url" ] && [ "$url" != null ]; do # the answer comes in pages; follow each next link
page=$(curl -s "$url")
jq -r '.results[] | [.id, ([.sampleset.samples[].ancestry_broad] | join("; ")), ([.performance_metrics.effect_sizes[] | "\(.name_short) \(.estimate)"] | join("; "))] | @tsv' <<<"$page"
url=$(jq -r '.next' <<<"$page")
done
Only additive scores are accepted (no dosage_*_weight columns, no is_dominant or is_recessive rows); step 25 stops on any other.
What the Script Does Internally¶
- Scoring files. Downloads the GRCh38-harmonised file of each score from the PGS Catalog FTP when it is not in
${GENOME_DIR}/prs_scores/yet, checked against the md5 the catalog publishes beside it. A file whose#HmPOS_buildheader is notGRCh38is refused. There is no fallback to the author-reported file, which is often GRCh37 or rsID-only and would score the wrong positions without any visible sign. - Scores as pgsc_calc reads them.
bin/collect_summary.py prs-formatwrites each file as a custom GRCh38 scoring file:chr_nameandchr_positionfrom the harmonisedhm_chrandhm_pos, the effect allele, the other allele (fromother_allele, elsehm_inferOtherAllelewhen it names one allele), the weight, and the catalog's trait as its label. Rows the catalog could not place on GRCh38 are dropped, and so are rows on chrY or the mitochondrial genome. Rows on chrX are kept when you give the sex and dropped when you do not (the step prints how many of each): pgsc_calc converts your genotypes with plink2, which stops on a chrX record when it has no sex. Without the sex, PGS000662 (prostate cancer) loses 8 of its 269 rows this way; see chrX and the sex. pgsc_calc then needs no network and no liftover. - Genotypes.
- With step 3's gVCF (the default since step 3 writes one): every score position, and with the panel every common panel SNV (the site list
setup.shwrites), is genotyped from the gVCF.bcftools convert --gvcf2vcfturns each reference block over a position into a 0/0 call with the reference base; a position with no coverage (./.) or outside every block stays missing. A 0/0 record gets as its ALT the position's first allele (a score's effect or other allele, the panel's ALT) that is not the reference, so pgsc_calc can match it. These genotypes are kept asprs/pgsc_calc/target.vcf.gz. - Without a gVCF (an older run): the variant-only VCF is cut to the same positions, so every site where you match the reference is missing (see below). When not one position is left, pgsc_calc is not started and every score is reported unmatched.
Either way pgsc_calc gets only those positions, a small file to convert: the autosomes, and chrX when you gave the sex.
4. pgsc_calc. With the sex, the step writes it for plink2 (pgsc_calc/sex/sex.tsv: the VCF's sample name, 1 for male, 2 for female) and hands it to pgsc_calc's PLINK2_VCF conversion as --update-sex, beside the arguments pgsc_calc gives that process itself. Runs pgsc_calc (from ${GENOME_DIR}/tools/pgsc_calc-<release>, which setup.sh or the step itself unpacks from GitHub's archive of the release, checked against PGSC_CALC_SHA256) with the images of versions.env, offline, its containers without network. With the panel it adds --run_ancestry. pgsc_calc matches each score's variants to your genotypes (strand flips, ambiguous A/T and C/G pairs dropped, one best match per variant), scores them with plink2 and, with the panel, projects you onto the panel's principal components and compares your score with the reference group most similar to you. A score that matches under 75% of its variants is dropped by pgsc_calc and gets no sum. With PRS_KEEP_AMBIGUOUS=true it adds --keep_ambiguous, and the ambiguous pairs are kept where step 2 left them in (see Strand-ambiguous rows).
5. Summary. bin/collect_summary.py prs-table reads pgsc_calc's match summary and scores into ${SAMPLE}_prs_summary.tsv, and with the panel writes step 26's ancestry table. pgsc_calc's work folder is deleted. Its results (its own HTML report, the match log) are kept, and so are its own trace, report and timeline: they move to pgsc_calc/pipeline_info/ as execution_trace.txt, execution_report.html and execution_timeline.html, without the time of the run in their names, and the rest of that folder (pgsc_calc's DAG) is deleted.
6. The previous run is kept. Before a rerun starts, the previous summary is copied to ${SAMPLE}_prs_summary.<YYYYMMDD-HHMM>.tsv (the time it was written) and the previous pgsc_calc/ folder is moved to pgsc_calc.<the same time>/, with -2, -3 and so on when that name is taken. Each run also copies its new summary to a dated name, with -2, -3 when an earlier summary of the same minute holds it. Two runs can then be compared side by side, and nothing is lost when a score update changes the baseline. Delete the old copies yourself when you no longer need them.
The Nextflow pipeline does the same with four processes: PRS_PREPARE, PRS_SCORE_SITES, PRS (pgsc_calc, which runs on the host because it starts its own containers) and PRS_SUMMARY. The sex is the samplesheet's sex column: a row with it is scored on chrX, a row without it on the autosomes, and PRS_PREPARE formats the scores once for each of the two that the run needs. Pass --ancestry_ref for the panel and --pgsc_calc ${GENOME_DIR}/tools/pgsc_calc-<release> to run offline; without it the PRS task fetches the release's archive, checked against PGSC_CALC_SHA256. PRS publishes pgsc_calc's output folder as prs/pgsc_calc/, with the same three files in its pipeline_info/. PRS_SUMMARY publishes the dated copy of the summary beside the summary, so a rerun keeps the earlier one; the published prs/pgsc_calc/ folder is replaced by the new run's, as Nextflow publishes over its own files. The bash step moves the old folder aside instead.
Which docker flags pgsc_calc's own containers get: they run as you (-u with your user and group ids) and with --network none, plus the sex folder's mount, in both entry points. In the pipeline they then get the run's own docker flags: its docker.runOptions and the PRS process's containerOptions (the process.containerOptions of the shared-host config, for example --cpu-shares 256), less their user (-u, --user) and network (--net, --network) flags, which stay the two above. The PRS task writes the result into images.config, the config it gives pgsc_calc, and prints it to its log. No other setting of the run reaches pgsc_calc's tasks: not its executor block, not a process setting, not another docker setting. pgsc_calc reads its own nextflow.config, images.config and Nextflow's user config ($NXF_HOME/config, by default ~/.nextflow/config), as any Nextflow run does, and images.config sets the docker.runOptions above. Step 25 runs outside any Nextflow run, so it has no outer flags to pass.
Output¶
| File | Contents |
|---|---|
${SAMPLE}_prs_summary.tsv |
Condition, PGS_ID, Score_SUM, Variants_Matched, Variants_Total, Matched_Pct, Percentile, Ancestry_Group, ChrX, Ambiguous, Applies_To, Sex_Check, Evaluated_In, Effect, Input. Evaluated_In and Effect are written when the score list names evaluations, as the shipped one does |
${SAMPLE}_prs_summary.<YYYYMMDD-HHMM>.tsv |
The same table under the time it was written, one per run, so a rerun keeps the earlier ones |
pgsc_calc.<YYYYMMDD-HHMM>/ |
Step 25 only: the previous run's pgsc_calc/ folder, moved aside by the rerun |
pgsc_calc/target.vcf.gz |
The genotypes pgsc_calc scored: the score (and panel) positions, from the gVCF or the VCF |
pgsc_calc/results/sample/score/ |
pgsc_calc's scores, its HTML report report.html, and with the panel the ancestry-adjusted scores and the principal components |
pgsc_calc/results/sample/match/ |
pgsc_calc's match log and summary |
pgsc_calc/pgsc_calc.log |
pgsc_calc's console output |
pgsc_calc/pipeline_info/ |
pgsc_calc's own Nextflow trace (execution_trace.txt: each of its tasks with its time, CPU and peak_rss), report and timeline. The trace holds task names, ids and figures only; the report and the timeline hold paths of your machine, as the pipeline's own pipeline_info/ does (Before you share outputs) |
../ancestry/${SAMPLE}_ancestry.tsv |
With the panel: step 26's table (population, its probability, principal components) |
All output is written to ${GENOME_DIR}/${SAMPLE}/prs/.
Runtime¶
Without the panel, 54 minutes to 1 h 17 min per genome in three measured 30-40x runs (runtime): the gVCF pass over every score position (the large scores reach most of the genome) took 21-31 minutes, and pgsc_calc 33-51 minutes with 24 GB. Preparing the scoring files took 6-11 minutes more, once per run. The panel adds its extraction, the intersection with your genotypes, the panel's PCA and the projection; see the measured numbers below.
Memory¶
Two tasks read nothing but the scoring files, so their memory follows the score list, not the genome, and a GitHub runner can measure it. On the nine scores the list had then (five of them have 1.7 to 6.9 million rows), tests/test_prs_memory.sh --measure (the E2E workflow's prs-memory job, on dispatch) measured on ubuntu-latest (4 CPUs, 16 GB RAM, containers limited to 14.1 GB):
| Task | Measured | Pipeline gives it | Step 25 gives it |
|---|---|---|---|
Scoring files written as pgsc_calc reads them (collect_summary.py prs-format, PRS_PREPARE) |
11.3 to 11.5 GB at the peak (memory.peak), about 5 minutes | 16 GB | no limit |
pgsc_calc's FORMAT_SCOREFILES (pgscatalog-format), one thread |
killed under the 14.1 GB limit without swap; with swap it finished, at about 14.2 GB of memory plus swap (sampled each second), in 17 to 21 minutes | one thread and all of PRS's 24 GB |
one thread and all of PGSC_MAX_MEMORY |
| the same, two threads (pgsc_calc's own choice) | about 15.7 GB, 13.5 minutes, with swap | ||
| the same, four threads (pgsc_calc's retry) | about 23.8 GB, 13 minutes, with swap |
In one whole pgsc_calc run on the PGS Catalog's synthetic genome-wide target (300 samples), with FORMAT_SCOREFILES on one thread and 15 GB, Nextflow's trace gave FORMAT_SCOREFILES a peak_rss of 18.1 GB (the sum of its processes' peaks, an upper bound), its two matching tasks 13.3 and 12.5 GB, and the run took 36 minutes. A real genome has one sample and fewer positions. On the three end-test genomes, with the current settings (PRS granted 24 GB), pgsc_calc passed on its first attempt each time: FORMAT_SCOREFILES held a working set of 12.7 to 13.5 GiB, the match combining 9 to 10 GiB and the variant matching 7 to 8.4 GiB, and PRS_PREPARE peaked at 11.3 to 12 GB (peak_rss).
A run of yours records the same figures for your genome and score list. Both entry points keep pgsc_calc's own trace as prs/pgsc_calc/pipeline_info/execution_trace.txt, with the peak_rss of each of its tasks (FORMAT_SCOREFILES, MATCH_VARIANTS, MATCH_COMBINE and the rest). The PRS row of the pipeline's own trace is only the host-side launcher, so read the nested trace to check or change the memory PRS asks for.
pgsc_calc v2.3.0 gives FORMAT_SCOREFILES 2 CPUs, and 4 on its own retry, under the same --max_memory, so its retry needs more memory, not less. Both entry points therefore give it one thread, the whole of pgsc_calc's memory, and no retry inside pgsc_calc at the same limit. One thread is the slowest of the three (4 to 8 minutes more) and the smallest.
- Pipeline:
PRSasks for 24 GB (conf/base.config), capped at--max_memory. With--max_memory 24.GBor more it gets all 24 GB. When a pgsc_calc task is killed for memory,PRSexits 137, so Nextflow retries it once with twice the memory, again capped at--max_memory.PRS_PREPAREasks for 16 GB, about 1.4 times its peak; 24 GB is about 1.3 times the 18.1 GB upper bound ofFORMAT_SCOREFILES. The labels gavePRS8 GB andPRS_PREPARE2 GB, and each was killed in a run on a ~30x sample (Lessons learned). - Step 25: pgsc_calc gets three quarters of the RAM, 24 GB on a 32 GB machine. When one of its tasks is killed for memory the step says so; run it again with a larger
PGSC_MAX_MEMORY. - A 16 GB machine caps both at about 12 to 15 GB, under what
FORMAT_SCOREFILESneeds. Docker lets a container swap as much again as its memory limit when the host has swap, which is how the whole pgsc_calc run above finished it with a 15 GB limit on a 16 GB runner; without swap it was killed.
Both figures depend on the score list: more or larger scoring files need more. The same job measures a changed list. The figures above were measured on the nine scores before the blood pressure and lipid scores were added. Those four have 101 to 884 rows each, against 1.7 to 6.9 million for the five large ones, so they add very little to the memory or the time of either task; they have not been measured on a whole genome yet. The asthma and height scores came later still. They add about 2.1 million rows to the roughly 24 million of the five large ones, so we expect both tasks to need about a tenth more memory. That is an estimate, not a measurement; the same job measures it.
The reference panel on a GitHub-hosted runner¶
The panel is pgsc_calc's 1000 Genomes database, pgsc_1000G_v1.tar.zst (PGSC_PANEL in versions.env), published by the PGS Catalog at https://ftp.ebi.ac.uk/pub/databases/spot/pgs/resources/. Measured by the E2E case tests/e2e/prs-3-panel-measure.sh on ubuntu-latest (4 CPUs, 16 GB RAM), with pgsc_calc projecting the PGS Catalog's synthetic genome-wide test target (600 samples) onto the panel and scoring PGS000018 (1.7 million variants):
| What | Measured (six dispatched runs) |
|---|---|
| Download | 7.43 GB (7,434,464,202 bytes). EBI served one connection at about 1.2 MB/s, about 100 minutes for the file; 8 byte ranges at once took 8 to 19 minutes |
| Site list | 61.6 million biallelic SNVs on the autosomes in the panel's GRCh38 table; setup.sh keeps the 6,966,553 with a panel frequency of 5% or more, the ones pgsc_calc's projection can use, in 2.7 minutes (plink2 one chromosome at a time, 8 GB; the whole table at once did not fit in 16 GB) |
| Disk while pgsc_calc runs | 14.9 to 16.1 GB on top of the panel (its work folder held 12.5 GB at the end): about 23 GB with the panel |
| Memory | 7.3 to 7.6 GB at the peak for the whole machine; the largest task, 6.9 to 7.1 GB (plink2 scoring the panel's own samples, and the panel's variant filter) |
| Time | 20 to 30 minutes for the run, on 4 CPUs |
That run scored PGS000018 alone and fits a 16 GB machine; the default scores need more (Memory). It does not fit the 14 GB of disk a GitHub-hosted runner guarantees, and the download alone outlasts the pull request test budget, so the panel stays opt-in (setup.sh --ancestry-panel), and the tests on every pull request project the fixture sample onto the PGS Catalog's 268 MB synthetic panel instead (tests/e2e/prs-2-ancestry.sh). The 1000 Genomes run above repeats on an E2E run dispatched on a branch other than main (or with PGSC_MEASURE=1): it adds about 40 minutes, which a dispatched run of main keeps for its other cases.
Interpreting Results¶
- Score_SUM: weighted sum of the effect alleles you carry. Higher = more genetic predisposition, but only relative to other people scored the same way.
- Variants_Matched: score variants pgsc_calc matched in your genotypes.
NAwhen every score fell under pgsc_calc's 75% match rate: pgsc_calc then stops and publishes no match counts, andMatched_Pctis the rate its log gave. - Variants_Total: variants of the scoring file with a GRCh38 position on the autosomes, plus chrX when the sex was given. Strand-ambiguous rows count unless a failed strand check left them out.
- Matched_Pct:
Variants_Matched / Variants_Totalas a percentage; the step warns below 50%, and pgsc_calc gives no sum below 75%. - Percentile: with the panel, where your score falls among the reference samples of
Ancestry_Group(pgsc_calc's empirical percentile,percentile_MostSimilarPop).NAwithout the panel, andNAfor a score that does not apply to the declared sex. - Ancestry_Group: the panel population whose genetic ancestry is most similar to yours (for the 1000 Genomes panel: AFR, AMR, EAS, EUR or SAS).
NAwithout the panel. - ChrX: the score's rows on chrX:
8 scoredwith the sex,8 left outwithout it,0when it has none. With the panel,8 left out (not matched in the panel run)when pgsc_calc matched none of them and5 scored, 3 left out (not matched in the panel run)when it matched some, see chrX and the sex. - Ambiguous: the score's strand-ambiguous rows (A/T and C/G) and what happened to them:
11 droppedby default,11 keptwithPRS_KEEP_AMBIGUOUS=true,11 dropped, strand check failedwhen the score's other rows did not agree with the reference.0when it has none. - Applies_To:
femaleormalefor a score built for one sex, empty for the rest. - Sex_Check: for those scores,
applies,not applicable to the declared sexorno sex declared; empty for the rest. - Evaluated_In: the groups the score has an evaluation for in the score list (
EUR;AFR;AMR), empty when it has none. - Effect: the OR or HR of the evaluation in
Ancestry_Group, with the group and the PPM id (HR 1.53 (EUR, PPM000591)). Empty without the panel, when the score has no percentile or was not evaluated in your group, and for a score that does not apply to the declared sex. - Input:
gvcfwhen the positions were genotyped from the gVCF,vcfwhen only the variant-only VCF was there.
Both reports show the percentile with its group, or one line saying "Raw score only" when no panel was used. With the panel they also show the effect in your group, or say that the score has none for it and name the groups it was evaluated in. When a score's chrX rows were left out because no sex was given, the same line names the score and how many rows; rows left out in a panel run show only in the summary's ChrX column. A score that does not apply to the declared sex shows "not applicable to the declared sex" in place of its percentile, and when no sex was declared the reports say which scores are for one sex only.
A raw Score_SUM from a run with the panel is not comparable with one from a run without it: with the panel pgsc_calc scores only the variants the panel has too. Compare percentiles, which only a panel run gives.
chrX and the sex¶
pgsc_calc converts your genotypes with plink2, and plink2 refuses chrX without the sample's sex. So a score's chrX rows are scored only when you give the sex (./scripts/25-prs.sh your_name male, or the samplesheet's sex column in Nextflow). Without it they are left out of the score: Variants_Total does not count them and the sum is the autosomes' alone. Of the 15 scores in the list, PGS000662 (prostate cancer) has 8 chrX rows of 269.
The two runs of the same score are therefore not comparable: with the sex the sum includes the chrX rows' weights, without it they are absent. Compare a sum only with sums made the same way. pgsc_calc's import does not split the pseudoautosomal regions out of chrX (no --split-par), so a row there is read with chrX's ploidy for the declared sex. The step does not take the sex from step 16's indexcov call: on a slice or a targeted BAM that call is not reliable (the test fixture reads a male sample as female), and a wrong sex would change the score silently.
With the ancestry panel, pgsc_calc keeps only the score variants that are also in the panel (its match step intersects the two, and the rest are marked excluded in the match log). A chrX row the panel lacks is then not scored even with the sex. Step 25 and PRS_SUMMARY read the match log in a panel run, and a score whose chrX rows were all excluded says N left out (not matched in the panel run) in its ChrX column, and one with some of them excluded names both counts (1 scored, 2 left out (not matched in the panel run)). The usual reason is that the panel lacks the site; a strand-ambiguous chrX row that pgsc_calc dropped reads the same, so the label names no cause. The e2e case tests/e2e/prs-scores-and-run-options-2-panel-chrx.sh checks this on the PGS Catalog's synthetic panel (GRCh38_HAPNEST_reference): pgsc_calc matched none of three chrX rows at sites the sample covers, so a man's prostate cancer score is scored on its autosomal rows only in a panel run. Its percentile is then taken among reference sums made the same way. The 1000 Genomes panel has not been checked this way yet.
Strand-ambiguous rows¶
A row whose two alleles are A and T, or C and G, reads the same on both strands, so a match cannot tell whether the study reported the forward strand or the reverse. pgsc_calc drops these rows by default, and so does step 25. Some scores lose many rows this way: of the 15, PGS002257 has 135 of 884.
PRS_KEEP_AMBIGUOUS=true (Nextflow --prs_keep_ambiguous) keeps them, under two conditions:
- WGS only. It needs step 3's gVCF, so every score position is genotyped against the reference. With only a variant VCF (or chip data) the step refuses it.
- A strand check per score.
prs-formatreads the reference base at each of the score's other SNV rows. The score passes when at least one row agrees with the reference on the forward strand and none needs a flip. A score that fails keeps its ambiguous rows out (Ambiguoussaysdropped, strand check failed), and the step prints the counts. PGS002257 fails it: one of its rows matches the reference only on the other strand, so its 135 ambiguous rows stay out even with the option. The three lipid scores pass.
The check says the score was reported on the forward strand; it does not prove each ambiguous row was. Keep it off when you want to compare with published sums, which mostly drop them too.
Hom-ref sites come from the gVCF¶
The VCF from step 3 lists only sites where you differ from the reference. A scoring variant whose effect allele is the reference allele is missing from that VCF when you carry two copies of it, and adds nothing to your sum. Step 3 also writes a gVCF, which records where you match the reference, and step 25 reads those sites from it as 0/0. With the gVCF, a site is missing only when it was not covered.
Without a gVCF (Input is vcf, a sample called by an older version), the sum misses the weight of every reference-allele effect allele you carry on two copies, and few scores reach pgsc_calc's 75% match rate. The step then prints:
NOTE: hom-ref sites are absent from this VCF, so the score is biased; not comparable to published distributions
Call the sample again with step 3 to get the gVCF. Imputation (step 14) is not needed for a score from a gVCF: it adds phasing, which the scores do not use, and a percentile needs the ancestry panel, not imputed sites.
What these scores are NOT¶
- A raw sum is NOT a percentile. Without the panel there is no percentile at all.
- A percentile is NOT a probability. Being at the 90th percentile does not mean a 90% chance of the condition.
- They are NOT comparable across conditions.
- The percentile compares you with the reference samples most similar to you, not with your own family or community. For a group the panel represents poorly, the percentile is less reliable; the ancestry table says when the population match is low-confidence.
- They are NOT stable across scoring file versions. A new harmonised file is a new baseline.
Limitations¶
- Most scores were developed in European-ancestry populations and predict less well for other ancestries, even with an ancestry-adjusted percentile.
- A score captures only the genetic component. Lifestyle, environment, and family history are often more predictive.
- Breast cancer and prostate cancer are scored for everyone, but get a percentile only for the sex they were built for, or when no sex was given.
- The test fixture is a few slices of one genome: its sums, match rates and percentiles show the code paths work, not what a real genome gives. Real percentiles, and the time and memory of the four new scores on a whole genome, have to come from a real run.
- Chip data matches far fewer score variants than a genome; most scores then fall under pgsc_calc's 75% match rate and get no sum (see Using chip data).
Notes¶
- Scoring files are cached in
${GENOME_DIR}/prs_scores/. Delete a file to download it again. A cached file without a#HmPOS_build=GRCh38header is deleted and downloaded again. - pgsc_calc starts with a fresh
pgsc_calc/folder each run: it keeps converted genotypes between runs, and an old run's would be scored instead of the current VCF. Step 25 moves the previous folder topgsc_calc.<YYYYMMDD-HHMM>/rather than deleting it.
Maintenance¶
- Recheck the PGS Catalog and pgsc_calc releases at least quarterly. A pgsc_calc bump means reading its
conf/modules.configat the new release and moving thePGSC_*_IMAGElines ofversions.envwith it, plusPGSC_CALC_SHA256andPGSC_CALC_NF_SCHEMA. - A scoring file update is a result-changing event: rerun step 25 and treat the output as a new baseline.
- Keep the
PGS ID, the harmonised file version, the pgsc_calc release and the pipeline commit together so score changes stay auditable.