Assess pLI in HPO genes

[bschilder]

Originally discussed here:

• Determining the appropriate background MSTExplorer#16 (comment)

We should asses the relationship between HPO genes and probability of Loss Intolerance (pLI) scores.

Questions

• Do HPO genes have greater pLI than randomly sampled background genes from the human genome?
• Is there a relationship between cell-type specificity scores (from our CTDs) and pLI?
• Do phenotypes with stronger cell type enrichment (in terms of fold-change) tend to have genes with greater pLI scores?
• Is there a relationship between pLI and GenCC evidence scores? (use HPOExplorer::get_gencc() to get the GenCC data)

@KittyMurphy it would be great if you could start exploring these questions and upload them as rmarkdown reports in this repo.
Perhaps start by just uploading a script showing me how to get the pLI scores for all genes, so I can start exploring this a bit as well.

Also just came across this paper:

• https://onlinelibrary.wiley.com/doi/10.1002/humu.23763

[bschilder]

https://onlinelibrary.wiley.com/doi/10.1002/humu.23763

Found this paper quite helpful in understanding pLI (and its many shortcomings as a metric).

Some additional (or alternative) metrics to consider:

• Variant Effect Predictor: now has a plugin for incorporating precomputed AlphaMissense scores! ensemblVEP is a bioc package for interfacing with VEP annotations.

[bschilder]

Ok, so while the data in VEP is super useful, the VEP software to access these annotations is a hot garbage fire that is virtually uninstallable. I tried everything under the sun:
https://gist.github.com/bschilder/8a64d266e0e3ab18075274ad539985ac

However, I was able to extract the AlphaMissense predictions directly! Turns out they already computed per gene scores for the entire protein-coding genome here (see their README):

AlphaMissense_gene_hg19.tsv.gz, AlphaMissense_gene_hg38.tsv.gz
Gene-level average predictions, which were computed by taking the mean
alphamissense_pathogenicity over all possible missense variants in a transcript
(canonical transcript).

With a little extra postprocessing, I got the gene symbols:

 am <- data.table::fread("https://storage.googleapis.com/dm_alphamissense/AlphaMissense_gene_hg38.tsv.gz")
  am$enst_id <- stringr::str_split(am$transcript_id,"\\.", simplify = TRUE)[,1] 
  map <- orthogene::map_genes(genes = unique(am$enst_id), 
                              target = "ENST",
                              species="human",
                              drop_na = FALSE, 
                              mthreshold = Inf)
  am_mapped <- unique(map[,c("input","name")]) |>
    data.table::data.table(key = "input") |>
    data.table::merge.data.table(am, by.x = "input", by.y = "enst_id")
am_mapped

pLI is still problably worth looking at, but I think ML-based AlphaMissense metric circumvents many of the shortcomings of the rule-based pLI metric.

[bschilder]

Found the latest pLI data from gnomad as well:
https://gnomad.broadinstitute.org/downloads/#v4-constraint

Importing that now for comparison with AlphaMissense.

readme <- suppressWarnings(
    readLines("https://storage.googleapis.com/gcp-public-data--gnomad/release/v4.0/constraint/README.txt")
  )
  pli <- data.table::fread("https://storage.googleapis.com/gcp-public-data--gnomad/release/v4.0/constraint/gnomad.v4.0.constraint_metrics.tsv")
  data.table::setorderv(pli, "mane_select",order=-1)

mane <- pli[mane_select==TRUE, lapply(.SD, mean, na.rm=TRUE), 
              .SDcols = is.numeric, by="gene"][, mane_select:=TRUE]
  pli_agg <- data.table::rbindlist(
    list(
      mane,
      pli[!gene %in% mane$gene, lapply(.SD, mean, na.rm=TRUE), 
          .SDcols = is.numeric, by="gene"][, mane_select:=FALSE]
    )
  )

[NathanSkene]

How would an AI model give us population frequency? I thought current generation AI models of protein folding are also really bad at predicting variant effects?

[KittyMurphy]

Just having a look into what I did previously when looking at pLI and genes under selective pressure.

I used the pLI for human transcripts from this study: The mutational constraint spectrum quantified from variation in 141,456 humans.

I'm attaching the relevant supplementary table:
supplementary_dataset_11_full_constraint_metrics.tsv.zip.

But as you've shared @bschilder, there is a more up to date version of pLI data.

[bschilder]

How would an AI model give us population frequency? I thought current generation AI models of protein folding are also really bad at predicting variant effects?

@NathanSkene Are you talking about variant population frequency, phenotype frequency, or disease frequency?
In any case, none of these were the intended usage of pLI as outlined here:
#50 (comment)

Also see here for my explanation of why pLI would not be appropriate for estimating population prevalence. Instead, getting epidemiological stats on population prevalence would make much more sense.

• Filter by frequency RareDiseasePrioritisation#34 (comment)