0.1 Data

0.1.1 UK Biobank Pharma Proteomics Project (UKB-PPP) Proteins

2940 summary statistics; Europeans; both GRCh19/38

Inflammation: 736
Cardiometabolic: 736
Oncology: 735
Neurology: 733

#!/bin/sh

# Directory containing the files
dir_path="/Users/charleenadams/ukbppp"

# Temporary file for storing category counts
temp_counts=$(mktemp)

# Loop through all .tar files in the directory
for file in "$dir_path"/*.tar; do
    # Extract the category from the filename
    category=$(echo "$file" | sed -E 's/.*_v1_([A-Za-z]+)(_II)?\.tar/\1/')

    # Output the category to the temp file
    echo "$category" >> "$temp_counts"
done

# Count the occurrences of each category and save to output
sort "$temp_counts" | uniq -c | awk '{print $1, $2}' | sort -nr > "$dir_path/category_counts.txt"

# Clean up
rm "$temp_counts"

echo "Category counts saved to $dir_path/category_counts.txt"

0.1.2 METSIM

1391 summary statistics; Finnish men; GRCh38; tsv files

Amino Acid: 215
Carbohydrate: 25
Cofactors and Vitamins: 38
Energy: 10
Lipid: 548
Nucleotide: 42
Partially Characterized: 16
Peptide: 42
Uncharacterized: 292
Xenobiotics: 163

0.1.3 Nightingale

249 summary statistics; Europeans; GRCh37; vcf files

Don’t have a breakdown by class, but do know the files include:

Lipoprotein Subclasses:
- Includes particle concentrations and composition
Lipoprotein Particle Size
Apolipoprotein A-I and B
Multiple Cholesterol and Triglyceride Measures
Albumin
Various Fatty Acids
Low-Molecular-Weight Metabolites:
- Amino acids (including branched-chain and aromatic)
- Glycolysis-related measures
- Ketone bodies

0.1.4 Clinical Phenotype(s)

Coronary Heart Disease (CHD) GRCh37; vcf file; unformatted and downloaded

I haven’t yet added CHD to the analysis, nor have I formatted the data. I’m waiting to finish benchmarking proteins and metabolites. But the plan is to add in CHD and other clinical phenotypes.

The CHD data below was done on the CARDIoGRAMplusC4D chip, which included Europeans and South Asians. (The authors of Hyprcoloc used it while saying that their method relies on the traits in the model being from the same population! Begs some questions in the colloquial sense [vs petitio principii], such as: “Which is it? Can Hyprcoloc handle traits from different ancestries or not?” “Was it on oversight or on purpose?”)

In the R documentation for Hyprcoloc, LD matrices aren’t required. This is because the authors assume the traits don’t overlap. There are options for indicating the traits do overlap, but the analysis gets vastly more complex. We need three additional matrices if samples overlap. So, it is best to choose wisely from the beginning, at least initally.

UKP-PPP, METSIM, and CHD CARDIoGRAMplusC4D should be reasonably independent. (I may need to change course on the UKB Nightingale idea. There are Nightingale summary statistics that contain fewer people that aren’t from the UKB. I might be able to get them instead of the ones I downloaded. METSIM is the priority, now.)

Phenotype	Abbreviation	Dataset	Author	PubMed ID	Sample Size (N)	Number of Cases (N Cases)	Population
Coronary Heart Disease	CHD	ieu-a-7	Nikpay	26343387	184,305	60,801	Mixed

NOTE: Remember the Jurgens data…

e.g, for BMI: “/Users/charleenadams/acy1_bmi/GWAS_sumstats_EUR__invnorm_bmi__TOTALsample.tsv”

0.2 Directories

Main: /Users/charleenadams/coloc_clin_mol
- Scripts: /Users/charleenadams/coloc_clin_mol/scripts
UKB-PPP files:
- Original .tar files were at /Users/charleenadams/ukbppp; deleted after extraction to save 1.5TB
- Extracted tars: /Users/charleenadams/ukbppp/coloc_clin_mol_extracted (main directory for UKB-PPP files now)
- rsIDs added files: /Users/charleenadams/merged_ukbppp
METSIM files: /Users/charleenadams/metsim_files
Nightingale files: /Users/charleenadams/met_d_files
Selected clinical phenotype files: /Users/charleenadams/clinical_files

0.3 Fetch metabolite summary statistics

0.4 Format UKB-PPP pQTLs

I had previously obtained these (see https://rpubs.com/YodaMendel/1243451) for an example of how to programmatically get files from UKB-PPP.

0.4.1 Untar all 2940


#!/bin/bash

# Define the directory containing your tar files
tar_directory="/Users/charleenadams/ukbppp"

# Define the output directory for extracted files
output_directory="/Users/charleenadams/ukbppp/coloc_clin_mol_extracted"

# Define the log file for tracking progress
log_file="/Users/charleenadams/ukbppp/extraction.log"

# Create the output directory if it doesn't exist
mkdir -p "$output_directory"

# Create or clear the log file
> "$log_file"

# Function to extract a single tar file
extract_tar() {
    local tar_file="$1"
    local output_directory="$2"
    local log_file="$3"

    # Log the start of extraction
    echo "Starting extraction: $tar_file" >> "$log_file"

    # Extract the tar file into the output directory
    if tar -xf "$tar_file" -C "$output_directory"; then
        echo "Completed extraction: $tar_file" >> "$log_file"
    else
        echo "Failed extraction: $tar_file" >> "$log_file"
    fi
}

export -f extract_tar

# Use `find` to locate all `.tar` files and process them in parallel
find "$tar_directory" -maxdepth 1 -type f -name "*.tar" | parallel -j 16 extract_tar {} "$output_directory" "$log_file"

echo "Extraction complete. Log saved to $log_file."

# Run with: 
#chmod +x untar.sh
#nohup ./untar.sh >untar.log 2>&1 &
ps -p 97839 -o stat,etime

0.4.2 Add rsIDs by CHR

~28 hours run time parallelized, including MAC sleeping overnight😞
- Next time use caffeinate -i ./add_rsids.sh

#!/bin/bash

# Base directory where your data is located
base_dir="/Users/charleenadams/ukbppp/coloc_clin_mol_extracted"

# Directory where RSID map files are located
map_dir="/Users/charleenadams/temp_BI/chemokine_rgs_olink/maps"

# Log file to track progress
log_file="$base_dir/add_rsids.log"
echo "Log started on $(date)" > "$log_file"

# Function to process each file
process_file() {
    local file="$1"
    local dir
    dir=$(dirname "$file")
    local chrom
    chrom=$(basename "$file" | sed -E 's/.*chr([0-9XY]+).*/\1/')
    local map_file="$map_dir/olink_rsid_map_mac5_info03_b0_7_chr${chrom}_patched_v2.tsv.gz"
    local output_file="${file%.gz}_rsid_added.gz"

    {
        echo "Processing file: $file (Chromosome $chrom)"

        if [[ -n "$chrom" && -f "$map_file" ]]; then
            echo "  Using RSID map: $map_file"

            # Remove existing rsid_added file to avoid overwrite prompts
            rm -f "$output_file"

            # Join the original file with the RSID map based on the ID column
            gunzip -c "$file" | \
            awk 'BEGIN {OFS="\t"}
                 NR==FNR {map[$1]=$4; genpos37[$1]=$5; genpos38[$1]=$6; next}
                 ($3 in map) {print $0, map[$3], genpos37[$3], genpos38[$3]}
                 !($3 in map) {print $0, ".", ".", "."}' <(gunzip -c "$map_file") - | \
            gzip -c > "$output_file"

            if [[ -s "$output_file" ]]; then
                echo "  RSID added with GENPOS37 and GENPOS38. Output saved to $output_file"

                # Delete the original file if the RSID-added file was successfully created
                rm -f "$file"
                echo "  Deleted original file: $file"
            else
                echo "  Warning: $output_file is empty! Original file not deleted."
            fi
        else
            echo "  Error: RSID map file for chromosome $chrom not found or chromosome not identified."
        fi
    } >> "$log_file" 2>&1
}

export -f process_file
export map_dir
export log_file

# Find all .gz files and process them in parallel
find "$base_dir" -type f -name "*.gz" | parallel -j16 process_file

echo "Log ended on $(date)" >> "$log_file"
echo "Processing complete. Check $log_file for details."

#Run with: 
#nano add_rsids.sh
#chmod +x add_rsids.sh
#nohup ./add_rsids.sh > add_rsids2.log 2>&1 &
#ps -p 54107 -o stat,etime

0.4.3 🚫 DON’T 🚫 merge UKP-PPP chromosome files: draft script!!!

Takes too long; keep as chromosomes; retained because it likely works and we might need the chromosomes merged if we try out a non-cis-region approach later.

#!/bin/bash

# Define directories
base_dir="/Users/charleenadams/ukbppp/coloc_clin_mol_extracted"
output_dir="/Users/charleenadams/merged_ukbppp"
log_file="${output_dir}/merge_rsid_files.log"

# Create output directory if it doesn't exist
mkdir -p "$output_dir"

# Start log file
echo "Processing started on $(date)" > "$log_file"

# Function to process each subdirectory
process_directory() {
    local dir="$1"

    # Extract the subdirectory name for naming the merged file
    subdir_name=$(basename "$dir")
    merged_file="${output_dir}/${subdir_name}_merged_rsid_added.gz"

    echo "Processing directory: $dir" | tee -a "$log_file"

    # Find all *_rsid_added.gz files in the current directory
    rsid_files=("${dir}"/*_rsid_added.gz)

    # Check if there are files to merge
    if [[ ${#rsid_files[@]} -gt 0 ]]; then
        # Extract header from the first file
        gunzip -c "${rsid_files[0]}" | head -n 1 > temp_header.txt

        # Process each file, excluding the header, and append to a temporary content file
        for file in "${rsid_files[@]}"; do
            gunzip -c "$file" | tail -n +2 >> temp_content.txt
        done

        # Concatenate header and content into the final merged file
        cat temp_header.txt temp_content.txt | gzip -c > "$merged_file"

        # Clean up temporary files
        rm temp_header.txt temp_content.txt

        # Confirm if the merged file was created successfully and isn’t empty
        if [[ -s "$merged_file" ]]; then
            echo "  Merged file created: $merged_file" | tee -a "$log_file"
        else
            echo "  Warning: Merged file is empty for directory $dir" | tee -a "$log_file"
        fi
    else
        echo "  No '_rsid_added.gz' files found in directory $dir" | tee -a "$log_file"
    fi
}

export -f process_directory
export log_file output_dir

# Find all subdirectories in base_dir and process them in parallel
find "$base_dir" -mindepth 1 -maxdepth 1 -type d -print0 | parallel -0 -j 15 process_directory

echo "Processing completed on $(date)" >> "$log_file"
echo "Check $log_file for details."

# Run with:
# Script located at: /Users/charleenadams/coloc_clin_mol/scripts
#nano chrom_merge.sh
#chmod +x chrom_merge.sh
#nohup ./chrom_merge.sh > chrom_merge.log 2>&1 &
#ps -p 20853 -o stat,etime

0.4.4 Create cis-regions

The method I devised below obtains cis-regions (500KB up and downstream of TSSes using Ensembl) for 2808 of 2940 (96%) of the files.

#!/bin/bash

# Directories
base_dir="/Users/charleenadams/ukbppp/coloc_clin_mol_extracted"
results_dir="/Users/charleenadams/cis_ukbppp"
mkdir -p "$results_dir"

# Log file
log_file="$results_dir/process_files.log"
echo "Processing started on $(date)" > "$log_file"

# Function to process each directory
process_directory() {
    local dir="$1"
    local protein_name
    protein_name=$(basename "$dir" | cut -d '_' -f1) # Extract protein name (e.g., ITGB6)

    echo "Processing directory: $dir (Protein: $protein_name)" | tee -a "$log_file"

    # Retrieve TSS for the protein using Ensembl REST API
    tss_info=$(curl -s "https://rest.ensembl.org/lookup/symbol/human/${protein_name}?content-type=application/json")
    if [[ -z "$tss_info" ]]; then
        echo "  Error: TSS not found for $protein_name" | tee -a "$log_file"
        return
    fi
    tss=$(echo "$tss_info" | jq '.start') # Extract TSS (start position)
    chrom=$(echo "$tss_info" | jq -r '.seq_region_name') # Extract chromosome

    # Check if TSS and chromosome are valid
    if [[ -z "$tss" || -z "$chrom" ]]; then
        echo "  Error: Invalid TSS or chromosome for $protein_name" | tee -a "$log_file"
        return
    fi

    # Define the region (500kb upstream and downstream)
    start=$((tss - 500000))
    end=$((tss + 500000))
    if ((start < 0)); then start=0; fi

    # Locate the chromosome file for the TSS
    chrom_file=$(find "$dir" -type f -name "discovery_chr${chrom}_*.gz")
    if [[ -z "$chrom_file" ]]; then
        echo "  Error: Chromosome file not found for $protein_name (chr $chrom)" | tee -a "$log_file"
        return
    fi

    # Filter the chromosome file for the region and save results
    output_file="${results_dir}/${protein_name}_cis_region.gz"
    gunzip -c "$chrom_file" | awk -v chrom="$chrom" -v start="$start" -v end="$end" \
        'NR==1 || ($1 == chrom && $2 >= start && $2 <= end)' | gzip -c > "$output_file"

    if [[ -s "$output_file" ]]; then
        echo "  Cis-region file created: $output_file" | tee -a "$log_file"
    else
        echo "  Warning: Cis-region file is empty for $protein_name" | tee -a "$log_file"
    fi
}

export -f process_directory
export results_dir log_file

# Run the process in parallel for all directories
find "$base_dir" -mindepth 1 -maxdepth 1 -type d | parallel -j15 process_directory {}

echo "Processing completed on $(date)" >> "$log_file"
echo "Check $log_file for details."

# Run with:
#nano cis_regions.sh
#chmod +x cis_regions.sh
#nohup ./cis_regions.sh > cis_regions.log 2>&1 &
#ps -p 20853 -o stat,etime

0.5 Debugging headers

Adding in the rsIDs created headers with mixed delimiters, so I fixed that.

#!/bin/bash

input_dir="/Users/charleenadams/cis_ukbppp"
output_dir="/Users/charleenadams/cis_ukbppp/cleaned"

mkdir -p "$output_dir"

for file in "$input_dir"/*_cis_region.gz; do
    if [[ -f "$file" ]]; then
        output_file="$output_dir/$(basename "$file" .gz)_cleaned.gz"
        echo "Processing $file -> $output_file"
        gunzip -c "$file" | \
        awk 'BEGIN {OFS="\t"} {for (i=1; i<=NF; i++) $i=$i} 1' | \
        gzip > "$output_file"
    else
        echo "Skipping $file: Not a valid file."
    fi
done

echo "Processing complete. Cleaned files are saved in $output_dir."

du -sh /Users/charleenadams/cis_ukbppp/* | grep -v '4.0K' | wc -l

du -sh /Users/charleenadams/cis_ukbppp/* | grep -v '4.0K' > non_4K_files.txt

0.6 Format Nightingale VCF files: draft script

#!/bin/bash

# Define the directory containing the .vcf.gz files
input_dir="/Users/charleenadams/met_d_files"
output_dir="/Users/charleenadams/met_d_files/formatted"
mkdir -p "$output_dir"

# Log file for tracking progress
log_file="$output_dir/format_vcf.log"
echo "Log started on $(date)" > "$log_file"

# Function to process each file
process_vcf() {
    local file="$1"
    local base_name
    base_name=$(basename "$file" .vcf.gz)
    local output_file="${output_dir}/${base_name}_formatted.vcf.gz"

    echo "Processing $file..." | tee -a "$log_file"

    # Extract the header and process the data
    gunzip -c "$file" | \
    awk -v OFS="\t" '
        BEGIN {process = 0}
        /^#CHROM/ {
            process = 1;
            print "CHROM", "POS19", "ID", "ALLELE0", "ALLELE1", "QUAL", "FILTER", "INFO", "BETA", "SE", "LOG10P", "A1FREQ", "rsid", "GENPOS";
        }
        process && !/^#/ {
            split($9, format_keys, ":");
            split($10, format_values, ":");
            # Create new columns based on FORMAT field
            for (i = 1; i <= length(format_keys); i++) {
                if (format_keys[i] == "ES") BETA = format_values[i];
                if (format_keys[i] == "SE") SE = format_values[i];
                if (format_keys[i] == "LP") LOG10P = format_values[i];
                if (format_keys[i] == "AF") A1FREQ = format_values[i];
                if (format_keys[i] == "ID") rsid = format_values[i];
            }
            # Map columns to required headers
            print $1, $2, rsid, $4, $5, $6, $7, $8, BETA, SE, LOG10P, A1FREQ, rsid, $2;
        }' | gzip -c > "$output_file"

    # Confirm processing success
    if [[ -s "$output_file" ]]; then
        echo "  Successfully processed: $output_file" | tee -a "$log_file"
        rm "$file"
    else
        echo "  Warning: Processing failed for $file" | tee -a "$log_file"
    fi
}

export -f process_vcf
export output_dir
export log_file

# Find all .vcf.gz files (excluding .tbi and .tbi.1 files) and process them in parallel
find "$input_dir" -type f -name "*.vcf.gz" ! -name "*.tbi*" | parallel -j4 process_vcf

echo "Log ended on $(date)" >> "$log_file"
echo "Processing complete. Check $log_file for details."

1 PCSK9 and METSIM

1.1 Subset METSIM files by PCSK9 cis-region

#!/bin/bash

#run with: 
#cd /Users/charleenadams/coloc_clin_mol/scripts
#chmod +x subset_PCSK9_metsim.sh
#nohup ./subset_PCSK9_metsim.sh > subset_PCSK9_metsim.log 2>&1 &
#ps -p 98320 -o stat,etime

# Input files and directories
cis_file="/Users/charleenadams/cis_ukbppp/cleaned/PCSK9_cis_region_cleaned.gz"
metsim_dir="/Users/charleenadams/metsim_files"
output_dir="/Users/charleenadams/metsim_files/metsim_cis_PCSK9"
log_file="subset_PCSK9_metsim.log"

# Create output directory
mkdir -p "$output_dir"

# Extract rsid list from the cis-region file
rsid_list=$(gzcat "$cis_file" | awk 'NR > 1 {print $1}' | sort | uniq)

if [[ -z "$rsid_list" ]]; then
  echo "Error: No rsid entries found in $cis_file" >> "$log_file"
  exit 1
fi

# Export variables for parallel
export rsid_list metsim_dir output_dir log_file

# Function to process a single METSIM file
process_metsim_file() {
  local metsim_file=$1
  local base_name=$(basename "$metsim_file" .txt.gz)
  local output_file="$output_dir/${base_name}_cis_PCSK9.txt"

  gzcat "$metsim_file" | awk -v rsids="$rsid_list" '
    BEGIN {
      split(rsids, rsid_array, " ");
      for (i in rsid_array) lookup[rsid_array[i]] = 1;
    }
    NR == 1 {print $0; next} # Print header
    $1 in lookup {print $0}
  ' > "$output_file"

  if [[ ! -s "$output_file" ]]; then
    echo "Warning: No matching rows for $metsim_file" >> "$log_file"
    rm "$output_file" # Remove empty files
  else
    echo "Processed $metsim_file" >> "$log_file"
  fi
}

export -f process_metsim_file

# Use GNU parallel to process files
find "$metsim_dir" -name "*.gz" | parallel -j $(nproc) process_metsim_file {}

1.2 Filter METSIM files on at least one SNP with P<5E-08 in PCSK9 cis-region

#!/usr/bin/env Rscript

# nohup Rscript ./filter_PCSK9_metsim_pval.R > filter_PCSK9_metsim_pval.log 2>&1 &
# ps -p 76539 -o stat,etime
# tail -f filter_PCSK9_metsim_pval.log

# Load necessary libraries
library(data.table)

# Directory where the txt files reside
directory <- "/Users/charleenadams/metsim_files/metsim_cis_PCSK9"
# Subdirectory where to copy files with pval < 5e-8
subdirectory <- "p_val_metsim_cis_PCSK9"

# Ensure subdirectory exists
if (!dir.exists(file.path(directory, subdirectory))) {
  dir.create(file.path(directory, subdirectory))
}

# Get all .txt files in the directory
files <- list.files(path = directory, pattern = "\\.txt$", full.names = TRUE)

# Process each file
for (file in files) {
  cat("Processing file:", basename(file), "\n")
  
  # Read the file
  data <- fread(file)
  
  # Check if any p-value is less than 5e-8
  if (any(data$pval < 5e-8, na.rm = TRUE)) {
    cat("File", basename(file), "has pval < 5e-8\n")
    file.copy(from = file, to = file.path(directory, subdirectory, basename(file)))
    cat("Copied", basename(file), "to", subdirectory, "\n")
  } else {
    cat("No match found in", basename(file), "\n")
  }
}

cat("Script completed\n")

1.3 Harmonizing PCSK9 and METSIM metabolites

#nohup Rscript ./harmonize_PCSK9_metsim_pval_subsets_parallel.R > harmonize_PCSK9_metsim_pval_subsets_parallel.log 2>&1 &
#ps -p 17652 -o stat,etime

# Load required libraries
library(data.table)
library(TwoSampleMR)

# Custom logging function
custom_log <- function(level, message) {
  msg <- paste(format(Sys.time(), "%Y-%m-%d %H:%M:%S"), level, message)
  cat(msg, "\n", file = "benchmark_harm_twosampleMR.log", append = TRUE)
  cat(msg, "\n")  # Also print to console
}

# Directories
cis_region_dir <- "/Users/charleenadams/cis_ukbppp/cleaned/"
metsim_dir <- "/Users/charleenadams/metsim_files/metsim_cis_PCSK9/p_val_metsim_cis_PCSK9"
output_dir <- "/Users/charleenadams/harmonized_results/"

# Function to read and format data for TwoSampleMR, including chromosome and position
read_and_format_data <- function(file_path, is_cis) {
  tryCatch({
    # Read file and convert to data.frame
    df <- fread(file_path, sep = "\t", header = TRUE)
    df <- as.data.frame(df)  # Ensure it is a data.frame
    
    # Format data for TwoSampleMR
    formatted_df <- format_data(df, 
                                type = ifelse(is_cis, "exposure", "outcome"), 
                                snp_col = ifelse(is_cis, "rsid", "rsids"), 
                                beta_col = ifelse(is_cis, "BETA", "beta"), 
                                se_col = ifelse(is_cis, "SE", "sebeta"), 
                                effect_allele_col = ifelse(is_cis, "ALLELE1", "alt"), 
                                other_allele_col = ifelse(is_cis, "ALLELE0", "ref"), 
                                eaf_col = ifelse(is_cis, "A1FREQ", "maf"), 
                                pval_col = ifelse(is_cis, "LOG10P", "pval"),
                                chr_col = ifelse(is_cis, "CHROM", "chrom"), 
                                pos_col = ifelse(is_cis, "POS38", "pos"),
                                phenotype_col = basename(file_path))
    
    # Set exposure or outcome names
    if (is_cis) {
      formatted_df$id.exposure <- basename(file_path)
    } else {
      formatted_df$id.outcome <- basename(file_path)
    }
    
    return(formatted_df)
  }, error = function(e) {
    custom_log("ERROR", paste("Error reading file", file_path, ":", e$message))
    return(NULL)
  })
}

# Main function using TwoSampleMR
main <- function() {
  # Create output directory if it doesn't exist
  dir.create(output_dir, showWarnings = FALSE, recursive = TRUE)
  
  # Select the cis-region file
  cis_file <- file.path(cis_region_dir, "PCSK9_cis_region_cleaned.gz")
  
  # Get the list of METSIM subset files
  metsim_files <- list.files(metsim_dir, pattern = "*.txt", full.names = TRUE)
  
  custom_log("INFO", paste("Processing cis-region file:", cis_file))
  custom_log("INFO", paste("Processing METSIM files:", paste(metsim_files, collapse = ", ")))
  
  # Read and format the cis-region file
  cis_df <- read_and_format_data(cis_file, is_cis = TRUE)
  if (is.null(cis_df)) {
    custom_log("ERROR", "Failed to read cis-region file.")
    return()
  }
  
  # Harmonize each METSIM file with the cis-region file
  harmonized_results <- lapply(metsim_files, function(metsim_file) {
    metsim_df <- read_and_format_data(metsim_file, is_cis = FALSE)
    if (!is.null(metsim_df)) {
      tryCatch({
        harmonized <- harmonise_data(exposure_dat = cis_df, outcome_dat = metsim_df, action = 2)
        custom_log("INFO", paste("Harmonized", metsim_file, "with", nrow(harmonized), "rows."))
        return(harmonized)
      }, error = function(e) {
        custom_log("ERROR", paste("Error harmonizing data for", metsim_file, ":", e$message))
        return(NULL)
      })
    }
    return(NULL)
  })
  
  # Combine harmonized results
  harmonized_results <- harmonized_results[!sapply(harmonized_results, is.null)]
  if (length(harmonized_results) > 0) {
    final_results <- rbindlist(harmonized_results)
    output_file <- file.path(output_dir, "harmonized_results_PCSK9.txt")
    fwrite(final_results, file = output_file, sep = "\t")
    custom_log("INFO", paste("Saved harmonized results to", output_file))
  } else {
    custom_log("WARNING", "No harmonized results were generated.")
  }
}

# Run the main function
main()

1.4 Prep for Hyprcoloc

Tricky. Do NOT under any circumstances select “update all” when prompted to install newer versions.

# Follow these installation steps exactly
remotes::install_version("Rcpp", version = "1.0.12", repos = "http://cran.us.r-project.org")
remotes::install_version("RcppEigen", version = "0.3.3.9.3", repos = "http://cran.us.r-project.org")

Sys.setenv(CC = "/opt/homebrew/opt/gcc/bin/gcc-14")
Sys.setenv(CXX = "/opt/homebrew/opt/gcc/bin/g++-14")

devtools::install_github("jrs95/hyprcoloc", 
                         build_opts = c("--no-resave-data", "--no-manual"), 
                         build_vignettes = TRUE)

# Load the necessary libraries
library(hyprcoloc)

# Regression coefficients from ten studies (Traits 1-5, 6-8 & 9-10 colocalize)
betas_t <- hyprcoloc::test.betas
ses_t <- hyprcoloc::test.ses

# Trait names and SNP IDs
traits_t <- paste0("T", 1:10)
rsid_t <- rownames(betas_t)
str(rsid_t)

# Colocalization analysis
test <- hyprcoloc(betas_t, ses_t, trait.names = traits_t, snp.id = rsid_t, snpscores = TRUE)

1.5 Hyprcoloc of PCSK9 and 184 METSIM metabolites

We conducted a colocalization analysis using HyPrColoc to explore shared genetic loci influencing PCSK9 and 184 METSIM metabolites.

# Run interactively in Rstudio
library(data.table)
library(hyprcoloc)
library(openxlsx)
library(DT)
library(knitr)

# Load phenotype mapping data
phenotype_map <- fread("/Users/charleenadams/metsim_files/phenotypes.tsv")
print(names(phenotype_map))  # Confirm structure of phenotype_map

##  [1] "chrom"         "pos"           "ref"           "alt"          
##  [5] "rsids"         "nearest_genes" "pval"          "num_peaks"    
##  [9] "phenocode"     "category"      "phenostring"

# Define file path to harmonized results
file_path <- "/Users/charleenadams/harmonized_results/harmonized_results_PCSK9.txt"

# Read the harmonized data
data <- fread(file_path)

# Remove all characters after and including the first `_`
data$id.outcome <- sub("_.*", "", data$id.outcome)
data$id.exposure <- sub("_.*", "", data$id.exposure)

# Ensure both data frames are data.tables
setDT(data)
setDT(phenotype_map)

# Merge phenostring from phenotype_map into data using id.outcome and phenocode
data <- merge(data, phenotype_map[, .(phenocode, phenostring)], 
              by.x = "id.outcome", by.y = "phenocode", all.x = TRUE)

# Rename the phenostring column to METSIM_name
setnames(data, "phenostring", "METSIM_name")

# Remove asterisks (*) from the METSIM_name column
data$METSIM_name <- gsub("\\*", "", data$METSIM_name)
cat("Data read successfully. Total rows:", nrow(data), "\n")

## Data read successfully. Total rows: 889665

# Ensure the data contains the required columns
required_columns <- c("SNP", "id.exposure", "METSIM_name", "beta.exposure", "beta.outcome", "se.exposure", "se.outcome")
if (!all(required_columns %in% colnames(data))) {
  stop("The harmonized data is missing required columns!")
}

# Extract unique traits and common SNPs
traits <- unique(c(data$id.exposure, data$METSIM_name))
common_snps <- unique(data$SNP)
cat("Number of traits:", length(traits), "Number of SNPs:", length(common_snps), "\n")

## Number of traits: 185 Number of SNPs: 4850

# Filter data to include only common SNPs
data <- data[SNP %in% common_snps]

# Reshape data to create a long table for matrix population
long_data <- rbindlist(list(
  data[, .(SNP, trait = id.exposure, beta = beta.exposure, se = se.exposure)],
  data[, .(SNP, trait = METSIM_name, beta = beta.outcome, se = se.outcome)]
), use.names = TRUE)

# Ensure no duplicate SNP-trait combinations
long_data <- unique(long_data, by = c("SNP", "trait"))
cat("Rows in reshaped long_data:", nrow(long_data), "\n")

## Rows in reshaped long_data: 894515

# Create matrices for betas and ses
betas_mat <- dcast(long_data, SNP ~ trait, value.var = "beta", fill = NA)
ses_mat <- dcast(long_data, SNP ~ trait, value.var = "se", fill = NA)
snps <- betas_mat$SNP

# Convert matrices to the appropriate format for hyprcoloc
rownames(betas_mat) <- betas_mat$SNP
rownames(ses_mat) <- ses_mat$SNP

betas_mat <- as.matrix(betas_mat[, -1, with = FALSE])
ses_mat <- as.matrix(ses_mat[, -1, with = FALSE])

# Restore the SNP column as rownames
rownames(betas_mat) <- snps
rownames(ses_mat) <- snps

# Confirm rownames are correctly set
head(rownames(betas_mat))

## [1] "rs1000364" "rs1002778" "rs1003323" "rs1003767" "rs1003768" "rs1003769"

cat("Dimensions of beta matrix:", dim(betas_mat), "\n")

## Dimensions of beta matrix: 4850 185

cat("Dimensions of SE matrix:", dim(ses_mat), "\n")

## Dimensions of SE matrix: 4850 185

# Check for rows with missing values
rows_with_na_betas <- apply(betas_mat, 1, function(row) any(is.na(row)))
rows_with_na_ses <- apply(ses_mat, 1, function(row) any(is.na(row)))

# Filter matrices to keep only rows with no missing values
betas_mat <- betas_mat[!rows_with_na_betas, , drop = FALSE]
ses_mat <- ses_mat[!rows_with_na_ses, , drop = FALSE]

cat("Dimensions of beta matrix:", dim(betas_mat), "\n")

## Dimensions of beta matrix: 4671 185

cat("Dimensions of SE matrix:", dim(ses_mat), "\n")

## Dimensions of SE matrix: 4671 185

# This noted-out solution assumes no effect of the missing SNP on the metabolite, which may be false. So, I removed all rows with any missing values instead, which would risk eliminating a causal SNP if in a row with some missing values. Perhaps six of one, half a dozen of the other. The same candidate causal SNP is found regardless.
# betas_mat_filtered[is.na(betas_mat_filtered)] <- 0
# ses_mat_filtered[is.na(ses_mat_filtered)] <- 1e6  # Large SE for low confidence

# Run hyprcoloc
cat("Running hyprcoloc...\n")

## Running hyprcoloc...

betas <- betas_mat
ses <- ses_mat
traits_names <- colnames(betas_mat)
snp_ids <- rownames(betas_mat)

# Uniform priors = FALSE
results_unipriorsFALSE <- hyprcoloc(betas, ses, trait.names = traits_names, snp.id = snp_ids, snpscores = TRUE,uniform.priors=FALSE)

# Static table with kable
kable(
  results_unipriorsFALSE$results,
  caption = "HyprColoc Results with Uniform Priors = FALSE",
  align = 'c'
)

HyprColoc Results with Uniform Priors = FALSE
iteration	traits	posterior_prob	regional_prob	candidate_snp	posterior_explained_by_snp	dropped_trait
1	1,2-dipalmitoyl-GPC (16:0/16:0), 1-(1-enyl-palmitoyl)-2-oleoyl-GPC (P-16:0/18:1), 1-(1-enyl-palmitoyl)-2-palmitoyl-GPC (P-16:0/16:0), 1-palmitoyl-2-stearoyl-GPC (16:0/18:0), 1-palmityl-2-palmitoyl-GPC (O-16:0/16:0), 1-stearoyl-2-oleoyl-GPI (18:0/18:1), 1-stearoyl-GPC (18:0), N-palmitoyl-sphingosine (d18:1/16:0), PCSK9, X - 23641, arachidoylcarnitine (C20), behenoyl sphingomyelin (d18:1/22:0), cerotoylcarnitine (C26), cholesterol, glycosyl ceramide (d18:1/20:0, d16:1/22:0), glycosyl ceramide (d18:2/24:1, d18:1/24:2), glycosyl-N-palmitoyl-sphingosine (d18:1/16:0), hydroxypalmitoyl sphingomyelin (d18:1/16:0(OH)), lactosyl-N-nervonoyl-sphingosine (d18:1/24:1), lignoceroyl sphingomyelin (d18:1/24:0), myristoyl dihydrosphingomyelin (d18:0/14:0), palmitoyl dihydrosphingomyelin (d18:0/16:0), palmitoyl sphingomyelin (d18:1/16:0), palmitoyl-sphingosine-phosphoethanolamine (d18:1/16:0), sphingomyelin (d17:1/16:0, d18:1/15:0, d16:1/17:0), sphingomyelin (d18:0/20:0, d16:0/22:0), sphingomyelin (d18:1/20:0, d16:1/22:0), sphingomyelin (d18:1/24:1, d18:2/24:0), sphingomyelin (d18:1/25:0, d19:0/24:1, d20:1/23:0, d19:1/24:0), sphingomyelin (d18:2/16:0, d18:1/16:1)	0.4593	0.6101	rs11591147	1	NA
2	(S)-3-hydroxybutyrylcarnitine, 3-hydroxydecanoate, 3-hydroxydecanoylcarnitine, 3-hydroxyhexanoate, 3-hydroxylaurate, 3-hydroxyoctanoate, X - 15469, X - 16397, X - 21353, X - 22519, cis-4-decenoate (10:1n6), dodecadienoate (12:2), tetradecadienoate (14:2)	0.3068	0.5345	rs1616691	1	NA
3	None	NA	0.0708	NA	NA	X - 15503
4	None	NA	0.0441	NA	NA	adipoylcarnitine (C6-DC)
5	None	NA	0.0290	NA	NA	dimethylarginine (SDMA + ADMA)
6	None	NA	0.0071	NA	NA	2,3-dihydroxy-5-methylthio-4-pentenoate (DMTPA)
7	None	NA	0.2606	NA	NA	decanoylcarnitine (C10)
8	None	NA	0.0050	NA	NA	kynurenine
9	None	NA	0.0205	NA	NA	octanoylcarnitine (C8)
10	None	NA	0.2715	NA	NA	dodecanedioate (C12-DC)
11	None	NA	0.0487	NA	NA	2-stearoyl-GPI (18:0)
12	None	NA	0.0395	NA	NA	tetradecanedioate (C14-DC)
13	None	NA	0.0207	NA	NA	N-acetylkynurenine (2)
14	None	NA	0.0036	NA	NA	2-palmitoyl-GPE (16:0)
15	None	NA	0.0036	NA	NA	X - 12101
16	None	NA	0.0059	NA	NA	1-stearoyl-2-arachidonoyl-GPI (18:0/20:4)
17	None	NA	0.0029	NA	NA	O-sulfo-L-tyrosine
18	None	NA	0.0039	NA	NA	hexadecanedioate (C16-DC)
19	None	NA	0.0334	NA	NA	X - 13431
20	None	NA	0.0124	NA	NA	1-palmitoyl-GPI (16:0)
21	None	NA	0.0075	NA	NA	3-methylglutarylcarnitine (2)
22	None	NA	0.0103	NA	NA	undecenoylcarnitine (C11:1)
23	None	NA	0.0062	NA	NA	alpha-hydroxyisovalerate
24	None	NA	0.0064	NA	NA	arabonate/xylonate
25	None	NA	0.0075	NA	NA	cis-4-decenoylcarnitine (C10:1)
26	None	NA	0.0033	NA	NA	N-palmitoyltaurine
27	None	NA	0.0375	NA	NA	ximenoylcarnitine (C26:1)
28	None	NA	0.0040	NA	NA	dodecenedioate (C12:1-DC)
29	None	NA	0.0044	NA	NA	X - 17676
30	None	NA	0.0183	NA	NA	uracil
31	None	NA	0.0138	NA	NA	N-palmitoylglycine
32	None	NA	0.0033	NA	NA	X - 18921
33	None	NA	0.0083	NA	NA	1-palmitoyl-2-oleoyl-GPI (16:0/18:1)
34	None	NA	0.0072	NA	NA	caprylate (8:0)
35	None	NA	0.0074	NA	NA	X - 12844
36	None	NA	0.0031	NA	NA	nervonoylcarnitine (C24:1)
37	None	NA	0.0074	NA	NA	nonanoylcarnitine (C9)
38	None	NA	0.0070	NA	NA	2-hydroxy-4-(methylthio)butanoic acid
39	None	NA	0.0038	NA	NA	glutamine conjugate of C6H10O2 (2)
40	None	NA	0.0035	NA	NA	quinolinate
41	None	NA	0.0040	NA	NA	1-palmitoyl-2-arachidonoyl-GPI (16:0/20:4)
42	None	NA	0.0034	NA	NA	3-hydroxysebacate
43	None	NA	0.0053	NA	NA	X - 12411
44	None	NA	0.0146	NA	NA	X - 16570
45	None	NA	0.0043	NA	NA	1-arachidonoyl-GPI (20:4)
46	None	NA	0.0090	NA	NA	N,N,N-trimethyl-5-aminovalerate
47	None	NA	0.0081	NA	NA	2-aminooctanoate
48	None	NA	0.0042	NA	NA	N-linoleoyltaurine
49	None	NA	0.0028	NA	NA	2-hydroxy-3-methylvalerate
50	None	NA	0.0041	NA	NA	glutarylcarnitine (C5-DC)
51	None	NA	0.0053	NA	NA	erythritol
52	None	NA	0.0035	NA	NA	X - 18922
53	None	NA	0.0029	NA	NA	kynurenate
54	None	NA	0.0027	NA	NA	lignoceroylcarnitine (C24)
55	None	NA	0.0040	NA	NA	N-oleoyltaurine
56	None	NA	0.0037	NA	NA	X - 14939
57	None	NA	0.0033	NA	NA	X - 11381
58	None	NA	0.0060	NA	NA	oleoyl-linoleoyl-glycerol (18:1/18:2) [2]
59	None	NA	0.0030	NA	NA	gamma-CEHC glucuronide
60	None	NA	0.0035	NA	NA	1-stearoyl-GPI (18:0)
61	None	NA	0.0036	NA	NA	hexanoylcarnitine (C6)
62	None	NA	0.0029	NA	NA	imidazole lactate
63	None	NA	0.0039	NA	NA	thyroxine
64	None	NA	0.0034	NA	NA	glutamine conjugate of C6H10O2 (1)
65	None	NA	0.0031	NA	NA	1-stearoyl-2-linoleoyl-GPI (18:0/18:2)
66	None	NA	0.0030	NA	NA	N-oleoylserine
67	None	NA	0.0032	NA	NA	1-oleoyl-GPS (18:1)
68	None	NA	0.0025	NA	NA	hexanoylglycine
69	None	NA	0.0028	NA	NA	X - 24462
70	None	NA	0.0599	NA	NA	glycine conjugate of C10H12O2
71	None	NA	0.0060	NA	NA	3-hydroxydodecanedioate
72	None	NA	0.0036	NA	NA	10-undecenoate (11:1n1)
73	None	NA	0.0029	NA	NA	choline phosphate
74	None	NA	0.0026	NA	NA	S-methylcysteine
75	None	NA	0.0053	NA	NA	phosphate
76	None	NA	0.0064	NA	NA	arabitol/xylitol
77	None	NA	0.0042	NA	NA	X - 21607
78	None	NA	0.0029	NA	NA	hexanoylglutamine
79	None	NA	0.0033	NA	NA	X - 18888
80	None	NA	0.0037	NA	NA	glycine conjugate of C10H14O2 (1)
81	None	NA	0.0030	NA	NA	phenyllactate (PLA)
82	None	NA	0.0025	NA	NA	1-palmitoyl-2-linoleoyl-GPI (16:0/18:2)
83	None	NA	0.0029	NA	NA	3-hydroxyadipate
84	None	NA	0.0025	NA	NA	ribulonate/xylulonate/lyxonate
85	None	NA	0.0032	NA	NA	4-guanidinobutanoate
86	None	NA	0.0037	NA	NA	3,7-dimethylurate
87	None	NA	0.0027	NA	NA	21-hydroxypregnenolone monosulfate (1)
88	None	NA	0.0035	NA	NA	1-palmitoleoyl-GPE (16:1)
89	None	NA	0.0040	NA	NA	1-(1-enyl-palmitoyl)-2-linoleoyl-GPE (P-16:0/18:2)
90	None	NA	0.0097	NA	NA	2-hydroxyoctanoate
91	None	NA	0.0027	NA	NA	flavin adenine dinucleotide (FAD)
92	None	NA	0.0025	NA	NA	caproate (6:0)
93	None	NA	0.0049	NA	NA	isoleucylglycine
94	None	NA	0.0025	NA	NA	N-acetylneuraminate
95	None	NA	0.0051	NA	NA	dehydroepiandrosterone sulfate (DHEA-S)
96	None	NA	0.0026	NA	NA	X - 24949
97	None	NA	0.0029	NA	NA	oleoyl ethanolamide
98	None	NA	0.0029	NA	NA	S-methylcysteine sulfoxide
99	None	NA	0.0025	NA	NA	pseudouridine
100	None	NA	0.0043	NA	NA	N-methylhydroxyproline
101	None	NA	0.0026	NA	NA	oleoyl-oleoyl-glycerol (18:1/18:1) [1]
102	None	NA	0.0031	NA	NA	hypotaurine
103	None	NA	0.0027	NA	NA	X - 23680
104	None	NA	0.0029	NA	NA	beta-hydroxyisovalerate
105	None	NA	0.0029	NA	NA	3-hydroxyhexanoylcarnitine (1)
106	None	NA	0.0025	NA	NA	4-allylphenol sulfate
107	None	NA	0.0029	NA	NA	X - 24518
108	None	NA	0.0026	NA	NA	5-hydroxyindoleacetate
109	None	NA	0.0029	NA	NA	X - 21355
110	None	NA	0.0027	NA	NA	tetradecadienedioate (C14:2-DC)
111	None	NA	0.0038	NA	NA	ribitol
112	None	NA	0.0026	NA	NA	N-stearoylserine
113	None	NA	0.0041	NA	NA	indoleacetylglutamine
114	None	NA	0.0024	NA	NA	glutamine conjugate of C7H12O2
115	None	NA	0.0028	NA	NA	X - 23639
116	None	NA	0.0028	NA	NA	indoleacetoylcarnitine
117	None	NA	0.0028	NA	NA	X - 17325
118	None	NA	0.0077	NA	NA	dimethyl sulfone
119	None	NA	0.0028	NA	NA	deoxycholate
120	None	NA	0.0026	NA	NA	X - 17438
121	None	NA	0.0024	NA	NA	ethylmalonate
122	None	NA	0.0025	NA	NA	linoleoyl ethanolamide
123	None	NA	0.0025	NA	NA	diacylglycerol (14:0/18:1, 16:0/16:1) [1]
124	None	NA	0.0025	NA	NA	2’-O-methylcytidine
125	None	NA	0.0035	NA	NA	serine
126	None	NA	0.0023	NA	NA	1-stearoyl-2-meadoyl-GPC (18:0/20:3n9)
127	None	NA	0.0034	NA	NA	cystathionine
128	None	NA	0.0033	NA	NA	urate
129	None	NA	0.0025	NA	NA	4-methylhexanoylglutamine
130	None	NA	0.0025	NA	NA	S-1-pyrroline-5-carboxylate
131	None	NA	0.0023	NA	NA	1-palmitoyl-GPC (16:0)
132	None	NA	0.0022	NA	NA	pimeloylcarnitine/3-methyladipoylcarnitine (C7-DC)
133	None	NA	0.0037	NA	NA	X - 18899
134	None	NA	0.0024	NA	NA	argininate
135	None	NA	0.0024	NA	NA	suberoylcarnitine (C8-DC)
136	None	NA	0.0023	NA	NA	mannitol/sorbitol
137	None	NA	0.0023	NA	NA	uridine
138	None	NA	0.0022	NA	NA	N-stearoyl-sphingadienine (d18:2/18:0)
139	None	NA	0.0025	NA	NA	X - 24337
140	None	NA	0.0022	NA	NA	X - 18901
141	None	NA	0.0022	NA	NA	21-hydroxypregnenolone disulfate
142	None	NA	0.0022	NA	NA	gamma-glutamylglutamate
143	None	NA	0.0021	NA	NA	glucose

# Uniform priors = TRUE
results_unipriorsTRUE <- hyprcoloc(betas, ses, trait.names = traits_names, snp.id = snp_ids, snpscores = TRUE, uniform.priors=TRUE)

kable(
  results_unipriorsTRUE$results,
  caption = "HyprColoc Results with Uniform Priors = TRUE",
  align = 'c'
)

HyprColoc Results with Uniform Priors = TRUE
iteration	traits	posterior_prob	regional_prob	candidate_snp	posterior_explained_by_snp	dropped_trait
1	1,2-dipalmitoyl-GPC (16:0/16:0), 1-(1-enyl-palmitoyl)-2-oleoyl-GPC (P-16:0/18:1), 1-(1-enyl-palmitoyl)-2-palmitoyl-GPC (P-16:0/16:0), 1-palmitoyl-2-stearoyl-GPC (16:0/18:0), 1-palmityl-2-palmitoyl-GPC (O-16:0/16:0), 1-stearoyl-2-oleoyl-GPI (18:0/18:1), 1-stearoyl-GPC (18:0), N-palmitoyl-sphingosine (d18:1/16:0), PCSK9, X - 23641, arachidoylcarnitine (C20), behenoyl sphingomyelin (d18:1/22:0), cerotoylcarnitine (C26), cholesterol, glycosyl ceramide (d18:1/20:0, d16:1/22:0), glycosyl ceramide (d18:2/24:1, d18:1/24:2), glycosyl-N-palmitoyl-sphingosine (d18:1/16:0), hydroxypalmitoyl sphingomyelin (d18:1/16:0(OH)), lactosyl-N-nervonoyl-sphingosine (d18:1/24:1), lignoceroyl sphingomyelin (d18:1/24:0), myristoyl dihydrosphingomyelin (d18:0/14:0), palmitoyl dihydrosphingomyelin (d18:0/16:0), palmitoyl sphingomyelin (d18:1/16:0), palmitoyl-sphingosine-phosphoethanolamine (d18:1/16:0), sphingomyelin (d17:1/16:0, d18:1/15:0, d16:1/17:0), sphingomyelin (d18:0/20:0, d16:0/22:0), sphingomyelin (d18:1/20:0, d16:1/22:0), sphingomyelin (d18:1/24:1, d18:2/24:0), sphingomyelin (d18:1/25:0, d19:0/24:1, d20:1/23:0, d19:1/24:0), sphingomyelin (d18:2/16:0, d18:1/16:1)	0.7712	0.7792	rs11591147	1.0000	NA
2	(S)-3-hydroxybutyrylcarnitine, 3-hydroxydecanoate, 3-hydroxydecanoylcarnitine, 3-hydroxyhexanoate, 3-hydroxylaurate, 3-hydroxyoctanoate, X - 15469, X - 16397, X - 18921, X - 21353, X - 22519, cis-4-decenoate (10:1n6), dodecadienoate (12:2), tetradecadienoate (14:2)	0.7638	0.7864	rs1616691	1.0000	NA
3	2,3-dihydroxy-5-methylthio-4-pentenoate (DMTPA), 2-palmitoyl-GPE (16:0), X - 15503, adipoylcarnitine (C6-DC), dimethylarginine (SDMA + ADMA), kynurenine, pseudouridine	0.7604	0.7791	rs141072298	0.8198	NA
4	N-acetylkynurenine (2), cis-4-decenoylcarnitine (C10:1), decanoylcarnitine (C10), hexanoylcarnitine (C6), octanoylcarnitine (C8)	0.8019	0.8254	rs148238018	0.9829	NA
5	None	NA	0.0106	NA	NA	dodecanedioate (C12-DC)
6	None	NA	0.0008	NA	NA	tetradecanedioate (C14-DC)
7	None	NA	0.0001	NA	NA	1-stearoyl-2-arachidonoyl-GPI (18:0/20:4)
8	None	NA	0.0001	NA	NA	X - 12101
9	None	NA	0.0000	NA	NA	O-sulfo-L-tyrosine
10	None	NA	0.0001	NA	NA	hexadecanedioate (C16-DC)
11	None	NA	0.0017	NA	NA	2-stearoyl-GPI (18:0)
12	None	NA	0.0002	NA	NA	1-palmitoyl-GPI (16:0)
13	None	NA	0.0001	NA	NA	alpha-hydroxyisovalerate
14	None	NA	0.0004	NA	NA	uracil
15	None	NA	0.0001	NA	NA	arabonate/xylonate
16	None	NA	0.0001	NA	NA	X - 12844
17	None	NA	0.0002	NA	NA	2-hydroxy-4-(methylthio)butanoic acid
18	None	NA	0.0001	NA	NA	1-palmitoyl-2-oleoyl-GPI (16:0/18:1)
19	None	NA	0.0001	NA	NA	caprylate (8:0)
20	None	NA	0.0001	NA	NA	glutamine conjugate of C6H10O2 (2)
21	None	NA	0.0001	NA	NA	1-palmitoyl-2-arachidonoyl-GPI (16:0/20:4)
22	None	NA	0.0001	NA	NA	N-palmitoyltaurine
23	None	NA	0.0001	NA	NA	flavin adenine dinucleotide (FAD)
24	None	NA	0.0001	NA	NA	3-hydroxysebacate
25	None	NA	0.0001	NA	NA	X - 12411
26	None	NA	0.0003	NA	NA	N-palmitoylglycine
27	None	NA	0.0002	NA	NA	oleoyl-linoleoyl-glycerol (18:1/18:2) [2]
28	None	NA	0.0004	NA	NA	ximenoylcarnitine (C26:1)
29	None	NA	0.0001	NA	NA	dodecenedioate (C12:1-DC)
30	None	NA	0.0001	NA	NA	1-arachidonoyl-GPI (20:4)
31	None	NA	0.0001	NA	NA	quinolinate
32	None	NA	0.0001	NA	NA	X - 17676
33	None	NA	0.0002	NA	NA	2-aminooctanoate
34	None	NA	0.0007	NA	NA	X - 13431
35	None	NA	0.0001	NA	NA	N-linoleoyltaurine
36	None	NA	0.0001	NA	NA	nervonoylcarnitine (C24:1)
37	None	NA	0.0000	NA	NA	2-hydroxy-3-methylvalerate
38	None	NA	0.0001	NA	NA	undecenoylcarnitine (C11:1)
39	None	NA	0.0000	NA	NA	imidazole lactate
40	None	NA	0.0002	NA	NA	N,N,N-trimethyl-5-aminovalerate
41	None	NA	0.0001	NA	NA	X - 18899
42	None	NA	0.0003	NA	NA	X - 16570
43	None	NA	0.0001	NA	NA	erythritol
44	None	NA	0.0001	NA	NA	1-stearoyl-GPI (18:0)
45	None	NA	0.0002	NA	NA	3-methylglutarylcarnitine (2)
46	None	NA	0.0001	NA	NA	X - 21607
47	None	NA	0.0001	NA	NA	10-undecenoate (11:1n1)
48	None	NA	0.0001	NA	NA	4-guanidinobutanoate
49	None	NA	0.0001	NA	NA	X - 18922
50	None	NA	0.0001	NA	NA	lignoceroylcarnitine (C24)
51	None	NA	0.0001	NA	NA	nonanoylcarnitine (C9)
52	None	NA	0.0001	NA	NA	hexanoylglutamine
53	None	NA	0.0001	NA	NA	N-oleoyltaurine
54	None	NA	0.0001	NA	NA	S-methylcysteine sulfoxide
55	None	NA	0.0001	NA	NA	X - 14939
56	None	NA	0.0001	NA	NA	2-hydroxyoctanoate
57	None	NA	0.0001	NA	NA	3,7-dimethylurate
58	None	NA	0.0001	NA	NA	1-oleoyl-GPS (18:1)
59	None	NA	0.0001	NA	NA	oleoyl ethanolamide
60	None	NA	0.0001	NA	NA	glutamine conjugate of C6H10O2 (1)
61	None	NA	0.0000	NA	NA	X - 24462
62	None	NA	0.0001	NA	NA	glutarylcarnitine (C5-DC)
63	None	NA	0.0001	NA	NA	hypotaurine
64	None	NA	0.0001	NA	NA	thyroxine
65	None	NA	0.0000	NA	NA	3-hydroxyadipate
66	None	NA	0.0001	NA	NA	isoleucylglycine
67	None	NA	0.0001	NA	NA	gamma-CEHC glucuronide
68	None	NA	0.0000	NA	NA	4-allylphenol sulfate
69	None	NA	0.0013	NA	NA	glycine conjugate of C10H12O2
70	None	NA	0.0001	NA	NA	glycine conjugate of C10H14O2 (1)
71	None	NA	0.0001	NA	NA	phenyllactate (PLA)
72	None	NA	0.0001	NA	NA	1-stearoyl-2-linoleoyl-GPI (18:0/18:2)
73	None	NA	0.0001	NA	NA	ribitol
74	None	NA	0.0001	NA	NA	X - 11381
75	None	NA	0.0001	NA	NA	3-hydroxydodecanedioate
76	None	NA	0.0000	NA	NA	5-hydroxyindoleacetate
77	None	NA	0.0001	NA	NA	indoleacetylglutamine
78	None	NA	0.0001	NA	NA	1-palmitoleoyl-GPE (16:1)
79	None	NA	0.0000	NA	NA	choline phosphate
80	None	NA	0.0000	NA	NA	N-oleoylserine
81	None	NA	0.0000	NA	NA	caproate (6:0)
82	None	NA	0.0001	NA	NA	N-methylhydroxyproline
83	None	NA	0.0000	NA	NA	X - 24949
84	None	NA	0.0001	NA	NA	dehydroepiandrosterone sulfate (DHEA-S)
85	None	NA	0.0001	NA	NA	X - 18888
86	None	NA	0.0000	NA	NA	ribulonate/xylulonate/lyxonate
87	None	NA	0.0000	NA	NA	X - 24518
88	None	NA	0.0001	NA	NA	beta-hydroxyisovalerate
89	None	NA	0.0000	NA	NA	3-hydroxyhexanoylcarnitine (1)
90	None	NA	0.0000	NA	NA	kynurenate
91	None	NA	0.0001	NA	NA	1-(1-enyl-palmitoyl)-2-linoleoyl-GPE (P-16:0/18:2)
92	None	NA	0.0001	NA	NA	dimethyl sulfone
93	None	NA	0.0000	NA	NA	arabitol/xylitol
94	None	NA	0.0000	NA	NA	oleoyl-oleoyl-glycerol (18:1/18:1) [1]
95	None	NA	0.0000	NA	NA	glutamine conjugate of C7H12O2
96	None	NA	0.0000	NA	NA	21-hydroxypregnenolone monosulfate (1)
97	None	NA	0.0000	NA	NA	1-palmitoyl-2-linoleoyl-GPI (16:0/18:2)
98	None	NA	0.0000	NA	NA	X - 23680
99	None	NA	0.0000	NA	NA	X - 23639
100	None	NA	0.0001	NA	NA	hexanoylglycine
101	None	NA	0.0001	NA	NA	serine
102	None	NA	0.0000	NA	NA	S-1-pyrroline-5-carboxylate
103	None	NA	0.0000	NA	NA	X - 21355
104	None	NA	0.0002	NA	NA	phosphate
105	None	NA	0.0000	NA	NA	indoleacetoylcarnitine
106	None	NA	0.0000	NA	NA	S-methylcysteine
107	None	NA	0.0000	NA	NA	tetradecadienedioate (C14:2-DC)
108	None	NA	0.0000	NA	NA	X - 17325
109	None	NA	0.0000	NA	NA	linoleoyl ethanolamide
110	None	NA	0.0000	NA	NA	ethylmalonate
111	None	NA	0.0000	NA	NA	1-palmitoyl-GPC (16:0)
112	None	NA	0.0000	NA	NA	1-stearoyl-2-meadoyl-GPC (18:0/20:3n9)
113	None	NA	0.0001	NA	NA	deoxycholate
114	None	NA	0.0000	NA	NA	2’-O-methylcytidine
115	None	NA	0.0001	NA	NA	cystathionine
116	None	NA	0.0000	NA	NA	urate
117	None	NA	0.0000	NA	NA	N-acetylneuraminate
118	None	NA	0.0000	NA	NA	N-stearoylserine
119	None	NA	0.0000	NA	NA	4-methylhexanoylglutamine
120	None	NA	0.0000	NA	NA	mannitol/sorbitol
121	None	NA	0.0000	NA	NA	X - 17438
122	None	NA	0.0000	NA	NA	pimeloylcarnitine/3-methyladipoylcarnitine (C7-DC)
123	None	NA	0.0000	NA	NA	diacylglycerol (14:0/18:1, 16:0/16:1) [1]
124	None	NA	0.0000	NA	NA	uridine
125	None	NA	0.0000	NA	NA	suberoylcarnitine (C8-DC)
126	None	NA	0.0000	NA	NA	N-stearoyl-sphingadienine (d18:2/18:0)
127	None	NA	0.0000	NA	NA	argininate
128	None	NA	0.0000	NA	NA	X - 24337
129	None	NA	0.0000	NA	NA	X - 18901
130	None	NA	0.0000	NA	NA	21-hydroxypregnenolone disulfate
131	None	NA	0.0000	NA	NA	gamma-glutamylglutamate
132	None	NA	0.0000	NA	NA	glucose

# Without assuming uniform priors and with reg.thresh and align.thresh
results_uniform_priors_FALSE_reg.thresh.95_align.thresh.95 <- hyprcoloc(betas, ses, trait.names = traits_names, snp.id = snp_ids, uniform.priors = FALSE,snpscores = TRUE,
  reg.thresh = 0.95, align.thresh = 0.95
)

kable(
  results_uniform_priors_FALSE_reg.thresh.95_align.thresh.95$results,
  caption = "HyprColoc Results with Uniform Priors = FALSE and regional and alignment thresholds = 0.95",
  align = 'c'
)

HyprColoc Results with Uniform Priors = FALSE and regional and alignment thresholds = 0.95
iteration	traits	posterior_prob	regional_prob	candidate_snp	posterior_explained_by_snp	dropped_trait
1	1,2-dipalmitoyl-GPC (16:0/16:0), 1-(1-enyl-palmitoyl)-2-oleoyl-GPC (P-16:0/18:1), 1-(1-enyl-palmitoyl)-2-palmitoyl-GPC (P-16:0/16:0), 1-palmitoyl-2-stearoyl-GPC (16:0/18:0), 1-palmityl-2-palmitoyl-GPC (O-16:0/16:0), N-palmitoyl-sphingosine (d18:1/16:0), PCSK9, behenoyl sphingomyelin (d18:1/22:0), cholesterol, glycosyl ceramide (d18:1/20:0, d16:1/22:0), glycosyl ceramide (d18:2/24:1, d18:1/24:2), glycosyl-N-palmitoyl-sphingosine (d18:1/16:0), hydroxypalmitoyl sphingomyelin (d18:1/16:0(OH)), lactosyl-N-nervonoyl-sphingosine (d18:1/24:1), lignoceroyl sphingomyelin (d18:1/24:0), myristoyl dihydrosphingomyelin (d18:0/14:0), palmitoyl dihydrosphingomyelin (d18:0/16:0), palmitoyl sphingomyelin (d18:1/16:0), palmitoyl-sphingosine-phosphoethanolamine (d18:1/16:0), sphingomyelin (d17:1/16:0, d18:1/15:0, d16:1/17:0), sphingomyelin (d18:0/20:0, d16:0/22:0), sphingomyelin (d18:1/20:0, d16:1/22:0), sphingomyelin (d18:1/24:1, d18:2/24:0), sphingomyelin (d18:2/16:0, d18:1/16:1)	0.9687	0.9808	rs11591147	1	NA
2	(S)-3-hydroxybutyrylcarnitine, 3-hydroxydecanoate, 3-hydroxydecanoylcarnitine, 3-hydroxyhexanoate, 3-hydroxylaurate, 3-hydroxyoctanoate, X - 15469, X - 21353, X - 22519	0.9832	0.9884	rs1616691	1	NA
3	None	NA	0.0708	NA	NA	X - 15503
4	None	NA	0.0441	NA	NA	adipoylcarnitine (C6-DC)
5	None	NA	0.0290	NA	NA	dimethylarginine (SDMA + ADMA)
6	None	NA	0.0071	NA	NA	2,3-dihydroxy-5-methylthio-4-pentenoate (DMTPA)
7	None	NA	0.2606	NA	NA	decanoylcarnitine (C10)
8	None	NA	0.0050	NA	NA	kynurenine
9	None	NA	0.0036	NA	NA	X - 12101
10	None	NA	0.0205	NA	NA	octanoylcarnitine (C8)
11	None	NA	0.7867	NA	NA	cis-4-decenoate (10:1n6)
12	None	NA	0.0044	NA	NA	2-stearoyl-GPI (18:0)
13	None	NA	0.0059	NA	NA	1-stearoyl-2-arachidonoyl-GPI (18:0/20:4)
14	None	NA	0.2715	NA	NA	dodecanedioate (C12-DC)
15	None	NA	0.0029	NA	NA	O-sulfo-L-tyrosine
16	None	NA	0.0395	NA	NA	tetradecanedioate (C14-DC)
17	None	NA	0.0036	NA	NA	2-palmitoyl-GPE (16:0)
18	None	NA	0.0124	NA	NA	1-palmitoyl-GPI (16:0)
19	None	NA	0.0139	NA	NA	N-acetylkynurenine (2)
20	None	NA	0.0375	NA	NA	ximenoylcarnitine (C26:1)
21	None	NA	0.0183	NA	NA	uracil
22	None	NA	0.0334	NA	NA	X - 13431
23	None	NA	0.0039	NA	NA	hexadecanedioate (C16-DC)
24	None	NA	0.0103	NA	NA	undecenoylcarnitine (C11:1)
25	None	NA	0.0033	NA	NA	N-palmitoyltaurine
26	None	NA	0.1594	NA	NA	dodecadienoate (12:2)
27	None	NA	0.0075	NA	NA	3-methylglutarylcarnitine (2)
28	None	NA	0.0062	NA	NA	alpha-hydroxyisovalerate
29	None	NA	0.0064	NA	NA	arabonate/xylonate
30	None	NA	0.0075	NA	NA	cis-4-decenoylcarnitine (C10:1)
31	None	NA	0.0138	NA	NA	N-palmitoylglycine
32	None	NA	0.0030	NA	NA	X - 16397
33	None	NA	0.0058	NA	NA	X - 18921
34	None	NA	0.0040	NA	NA	1-palmitoyl-2-arachidonoyl-GPI (16:0/20:4)
35	None	NA	0.0072	NA	NA	caprylate (8:0)
36	None	NA	0.0074	NA	NA	X - 12844
37	None	NA	0.0040	NA	NA	dodecenedioate (C12:1-DC)
38	None	NA	0.0036	NA	NA	X - 17676
39	None	NA	0.0088	NA	NA	1-arachidonoyl-GPI (20:4)
40	None	NA	0.0070	NA	NA	2-hydroxy-4-(methylthio)butanoic acid
41	None	NA	0.0083	NA	NA	1-palmitoyl-2-oleoyl-GPI (16:0/18:1)
42	None	NA	0.0053	NA	NA	X - 12411
43	None	NA	0.0081	NA	NA	2-aminooctanoate
44	None	NA	0.0035	NA	NA	quinolinate
45	None	NA	0.0037	NA	NA	N-linoleoyltaurine
46	None	NA	0.0031	NA	NA	nervonoylcarnitine (C24:1)
47	None	NA	0.0090	NA	NA	N,N,N-trimethyl-5-aminovalerate
48	None	NA	0.0034	NA	NA	3-hydroxysebacate
49	None	NA	0.0146	NA	NA	X - 16570
50	None	NA	0.0074	NA	NA	nonanoylcarnitine (C9)
51	None	NA	0.0035	NA	NA	X - 18922
52	None	NA	0.0028	NA	NA	2-hydroxy-3-methylvalerate
53	None	NA	0.0029	NA	NA	kynurenate
54	None	NA	0.0060	NA	NA	oleoyl-linoleoyl-glycerol (18:1/18:2) [2]
55	None	NA	0.0030	NA	NA	gamma-CEHC glucuronide
56	None	NA	0.0035	NA	NA	1-stearoyl-GPI (18:0)
57	None	NA	0.0037	NA	NA	glutamine conjugate of C6H10O2 (2)
58	None	NA	0.0040	NA	NA	N-oleoyltaurine
59	None	NA	0.0041	NA	NA	glutarylcarnitine (C5-DC)
60	None	NA	0.0026	NA	NA	S-methylcysteine
61	None	NA	0.0033	NA	NA	X - 11381
62	None	NA	0.0036	NA	NA	hexanoylcarnitine (C6)
63	None	NA	0.0032	NA	NA	X - 18899
64	None	NA	0.0031	NA	NA	lignoceroylcarnitine (C24)
65	None	NA	0.0034	NA	NA	erythritol
66	None	NA	0.0034	NA	NA	glutamine conjugate of C6H10O2 (1)
67	None	NA	0.0032	NA	NA	1-oleoyl-GPS (18:1)
68	None	NA	0.0028	NA	NA	X - 24462
69	None	NA	0.0115	NA	NA	X - 21607
70	None	NA	0.0035	NA	NA	tetradecadienoate (14:2)
71	None	NA	0.0032	NA	NA	hexanoylglutamine
72	None	NA	0.0029	NA	NA	imidazole lactate
73	None	NA	0.0039	NA	NA	thyroxine
74	None	NA	0.0028	NA	NA	hypotaurine
75	None	NA	0.0131	NA	NA	arachidoylcarnitine (C20)
76	None	NA	0.0033	NA	NA	1-palmitoleoyl-GPE (16:1)
77	None	NA	0.0044	NA	NA	ribitol
78	None	NA	0.0033	NA	NA	4-guanidinobutanoate
79	None	NA	0.0037	NA	NA	X - 14939
80	None	NA	0.0043	NA	NA	N-methylhydroxyproline
81	None	NA	0.0060	NA	NA	3-hydroxydodecanedioate
82	None	NA	0.0029	NA	NA	3-hydroxyadipate
83	None	NA	0.0599	NA	NA	glycine conjugate of C10H12O2
84	None	NA	0.0033	NA	NA	X - 18888
85	None	NA	0.0030	NA	NA	N-oleoylserine
86	None	NA	0.0164	NA	NA	1-stearoyl-GPC (18:0)
87	None	NA	0.0026	NA	NA	caproate (6:0)
88	None	NA	0.0123	NA	NA	phosphate
89	None	NA	0.0049	NA	NA	isoleucylglycine
90	None	NA	0.0037	NA	NA	10-undecenoate (11:1n1)
91	None	NA	0.0031	NA	NA	sphingomyelin (d18:1/25:0, d19:0/24:1, d20:1/23:0, d19:1/24:0)
92	None	NA	0.0052	NA	NA	1-stearoyl-2-oleoyl-GPI (18:0/18:1)
93	None	NA	0.0025	NA	NA	hexanoylglycine
94	None	NA	0.0043	NA	NA	glycine conjugate of C10H14O2 (1)
95	None	NA	0.0030	NA	NA	phenyllactate (PLA)
96	None	NA	0.0026	NA	NA	oleoyl-oleoyl-glycerol (18:1/18:1) [1]
97	None	NA	0.0035	NA	NA	1-stearoyl-2-linoleoyl-GPI (18:0/18:2)
98	None	NA	0.0029	NA	NA	3-hydroxyhexanoylcarnitine (1)
99	None	NA	0.0035	NA	NA	5-hydroxyindoleacetate
100	None	NA	0.0037	NA	NA	3,7-dimethylurate
101	None	NA	0.0031	NA	NA	X - 24518
102	None	NA	0.0040	NA	NA	1-(1-enyl-palmitoyl)-2-linoleoyl-GPE (P-16:0/18:2)
103	None	NA	0.0027	NA	NA	X - 23641
104	None	NA	0.0097	NA	NA	2-hydroxyoctanoate
105	None	NA	0.0025	NA	NA	choline phosphate
106	None	NA	0.0025	NA	NA	1-palmitoyl-2-linoleoyl-GPI (16:0/18:2)
107	None	NA	0.0025	NA	NA	N-acetylneuraminate
108	None	NA	0.0027	NA	NA	X - 23680
109	None	NA	0.0025	NA	NA	ribulonate/xylulonate/lyxonate
110	None	NA	0.0034	NA	NA	oleoyl ethanolamide
111	None	NA	0.0026	NA	NA	X - 24949
112	None	NA	0.0051	NA	NA	dehydroepiandrosterone sulfate (DHEA-S)
113	None	NA	0.0028	NA	NA	X - 17438
114	None	NA	0.0029	NA	NA	S-methylcysteine sulfoxide
115	None	NA	0.0028	NA	NA	beta-hydroxyisovalerate
116	None	NA	0.0025	NA	NA	pseudouridine
117	None	NA	0.0028	NA	NA	X - 23639
118	None	NA	0.0028	NA	NA	X - 17325
119	None	NA	0.0030	NA	NA	cerotoylcarnitine (C26)
120	None	NA	0.0037	NA	NA	dimethyl sulfone
121	None	NA	0.0028	NA	NA	indoleacetoylcarnitine
122	None	NA	0.0028	NA	NA	arabitol/xylitol
123	None	NA	0.0033	NA	NA	serine
124	None	NA	0.0031	NA	NA	flavin adenine dinucleotide (FAD)
125	None	NA	0.0026	NA	NA	21-hydroxypregnenolone monosulfate (1)
126	None	NA	0.0025	NA	NA	S-1-pyrroline-5-carboxylate
127	None	NA	0.0027	NA	NA	tetradecadienedioate (C14:2-DC)
128	None	NA	0.0025	NA	NA	glutamine conjugate of C7H12O2
129	None	NA	0.0025	NA	NA	linoleoyl ethanolamide
130	None	NA	0.0024	NA	NA	ethylmalonate
131	None	NA	0.0031	NA	NA	deoxycholate
132	None	NA	0.0025	NA	NA	diacylglycerol (14:0/18:1, 16:0/16:1) [1]
133	None	NA	0.0024	NA	NA	1-palmitoyl-GPC (16:0)
134	None	NA	0.0026	NA	NA	2’-O-methylcytidine
135	None	NA	0.0033	NA	NA	urate
136	None	NA	0.0034	NA	NA	cystathionine
137	None	NA	0.0023	NA	NA	1-stearoyl-2-meadoyl-GPC (18:0/20:3n9)
138	None	NA	0.0028	NA	NA	indoleacetylglutamine
139	None	NA	0.0028	NA	NA	X - 21355
140	None	NA	0.0027	NA	NA	N-stearoylserine
141	None	NA	0.0025	NA	NA	4-methylhexanoylglutamine
142	None	NA	0.0025	NA	NA	4-allylphenol sulfate
143	None	NA	0.0022	NA	NA	pimeloylcarnitine/3-methyladipoylcarnitine (C7-DC)
144	None	NA	0.0024	NA	NA	argininate
145	None	NA	0.0024	NA	NA	suberoylcarnitine (C8-DC)
146	None	NA	0.0023	NA	NA	mannitol/sorbitol
147	None	NA	0.0023	NA	NA	uridine
148	None	NA	0.0022	NA	NA	N-stearoyl-sphingadienine (d18:2/18:0)
149	None	NA	0.0025	NA	NA	X - 24337
150	None	NA	0.0022	NA	NA	X - 18901
151	None	NA	0.0022	NA	NA	21-hydroxypregnenolone disulfate
152	None	NA	0.0022	NA	NA	gamma-glutamylglutamate
153	None	NA	0.0021	NA	NA	glucose

# Save the results as tabs in Excel
# Define file path for the output Excel workbook
output_file <- "/Users/charleenadams/pcsk9_hyprcoloc_corrected_results_summary_with_snpscores.xlsx"

# Initialize a new workbook
wb <- createWorkbook()

# Function to add a result set to the workbook
add_results_to_workbook <- function(workbook, result_data, sheet_name) {
  # Extract the results data frame
  results_df <- result_data$results
  
  # Create a worksheet and write the data
  addWorksheet(workbook, sheet_name)
  writeData(workbook, sheet_name, results_df)
  
  # Apply styling
  header_style <- createStyle(fontSize = 12, fontColour = "#FFFFFF", halign = "CENTER", 
                              fgFill = "#4F81BD", textDecoration = "Bold")
  addStyle(workbook, sheet_name, header_style, rows = 1, cols = 1:ncol(results_df), gridExpand = TRUE)
}

# Function to add SNP scores to the workbook with rsIDs
add_snpscores_to_workbook <- function(workbook, snp_scores, rsIDs, sheet_name) {
  # Convert SNP scores to a data table and include rsIDs
  snp_scores_df <- as.data.table(snp_scores)
  colnames(snp_scores_df) <- paste0("SNP_Score_Set_", seq_len(ncol(snp_scores_df)))
  snp_scores_df$rsID <- rsIDs
  
  # Reorder columns to put rsID first
  setcolorder(snp_scores_df, c("rsID", setdiff(names(snp_scores_df), "rsID")))
  
  # Create a worksheet and write the data
  addWorksheet(workbook, sheet_name)
  writeData(workbook, sheet_name, snp_scores_df)
  
  # Apply styling
  header_style <- createStyle(fontSize = 12, fontColour = "#FFFFFF", halign = "CENTER", 
                              fgFill = "#4F81BD", textDecoration = "Bold")
  addStyle(workbook, sheet_name, header_style, rows = 1, cols = 1:ncol(snp_scores_df), gridExpand = TRUE)
}

# Add results and SNP scores for each analysis
# Ensure rsIDs are available as rownames in the matrices
rsIDs <- rownames(betas_mat)

# results_unipriorsFALSE
add_results_to_workbook(wb, results_unipriorsFALSE, "Results_Unipriors_FALSE")
add_snpscores_to_workbook(wb, results_unipriorsFALSE$snpscores, rsIDs, "SNP_Scores_Unipriors_FALSE")

# results_unipriorsTRUE
add_results_to_workbook(wb, results_unipriorsTRUE, "Results_Unipriors_TRUE")
add_snpscores_to_workbook(wb, results_unipriorsTRUE$snpscores, rsIDs, "SNP_Scores_Unipriors_TRUE")

# results_uniform_priors_FALSE_reg.thresh.95_align.thresh.95
add_results_to_workbook(wb, results_uniform_priors_FALSE_reg.thresh.95_align.thresh.95, 
                        "Unipriors_FALSE_Thresh_95")
add_snpscores_to_workbook(wb, results_uniform_priors_FALSE_reg.thresh.95_align.thresh.95$snpscores, 
                          rsIDs, "SNP_Scores_Thresh_95")

# Save the workbook
saveWorkbook(wb, output_file, overwrite = TRUE)

# Confirm completion
cat("Results and SNP scores saved to:", output_file, "\n")

## Results and SNP scores saved to: /Users/charleenadams/pcsk9_hyprcoloc_corrected_results_summary_with_snpscores.xlsx

1.6 Hyprcoloc of PCSK9 and top 10 METSIM metabolites

top = 10 metabolites containing rsIDs with smallest p-values

# Run interactively in Rstudio
library(data.table)
library(hyprcoloc)
library(openxlsx)
library(DT)
library(knitr)

# Load phenotype mapping data
phenotype_map <- fread("/Users/charleenadams/metsim_files/phenotypes.tsv")
print(names(phenotype_map))  # Confirm structure of phenotype_map

##  [1] "chrom"         "pos"           "ref"           "alt"          
##  [5] "rsids"         "nearest_genes" "pval"          "num_peaks"    
##  [9] "phenocode"     "category"      "phenostring"

# Define file path to harmonized results
file_path <- "/Users/charleenadams/harmonized_results/harmonized_results_PCSK9.txt"

# Read the harmonized data
data <- fread(file_path)

# Remove all characters after and including the first `_`
data$id.outcome <- sub("_.*", "", data$id.outcome)
data$id.exposure <- sub("_.*", "", data$id.exposure)

# Ensure both data frames are data.tables
setDT(data)
setDT(phenotype_map)

# Merge phenostring from phenotype_map into data using id.outcome and phenocode
data <- merge(data, phenotype_map[, .(phenocode, phenostring)], 
              by.x = "id.outcome", by.y = "phenocode", all.x = TRUE)

# Rename the phenostring column to METSIM_name
setnames(data, "phenostring", "METSIM_name")

# Remove asterisks (*) from the METSIM_name column
data$METSIM_name <- gsub("\\*", "", data$METSIM_name)
cat("Data read successfully. Total rows:", nrow(data), "\n")

## Data read successfully. Total rows: 889665

# Ensure the data table is sorted by pval.outcome in ascending order
data <- data[order(pval.outcome)]

# Get the top 10 unique id.outcome with the smallest pval.outcome
top_id_outcomes <- unique(data$id.outcome[order(data$pval.outcome)])[1:10]

# Subset the data to include only rows with these id.outcome values
subset_data <- data[id.outcome %in% top_id_outcomes]

# Verify the subset
cat("Number of rows in the subset:", nrow(subset_data), "\n")

## Number of rows in the subset: 48384

cat("Unique id.outcome in the subset:", length(unique(subset_data$id.outcome)), "\n")

## Unique id.outcome in the subset: 10

cat("Metabolites with the smallest p-values:\n", paste(unique(subset_data$METSIM_name), collapse = "\n"), "\n")

## Metabolites with the smallest p-values:
##  palmitoyl dihydrosphingomyelin (d18:0/16:0)
## palmitoyl sphingomyelin (d18:1/16:0)
## palmitoyl-sphingosine-phosphoethanolamine (d18:1/16:0)
## hydroxypalmitoyl sphingomyelin (d18:1/16:0(OH))
## cholesterol
## sphingomyelin (d18:1/24:1, d18:2/24:0)
## glycosyl ceramide (d18:1/20:0, d16:1/22:0)
## X - 15469
## 1-(1-enyl-palmitoyl)-2-palmitoyl-GPC (P-16:0/16:0)
## glycosyl-N-palmitoyl-sphingosine (d18:1/16:0)

# Ensure the subset_data contains the required columns
required_columns <- c("SNP", "id.exposure", "METSIM_name", "beta.exposure", "beta.outcome", "se.exposure", "se.outcome")
if (!all(required_columns %in% colnames(subset_data))) {
  stop("The harmonized subset_data is missing required columns!")
}

# Extract unique traits and common SNPs
traits <- unique(c(subset_data$id.exposure, subset_data$METSIM_name))
common_snps <- unique(subset_data$SNP)
cat("Number of traits:", length(traits), "Number of SNPs:", length(common_snps), "\n")

## Number of traits: 11 Number of SNPs: 4850

# Filter subset_data to include only common SNPs
filtered_subset_data <- subset_data[SNP %in% common_snps]

# Reshape subset_data to create a long table for matrix population
long_subset_data <- rbindlist(list(
  subset_data[, .(SNP, trait = id.exposure, beta = beta.exposure, se = se.exposure)],
  subset_data[, .(SNP, trait = METSIM_name, beta = beta.outcome, se = se.outcome)]
), use.names = TRUE)

# Ensure no duplicate SNP-trait combinations
long_subset_data <- unique(long_subset_data, by = c("SNP", "trait"))
cat("Rows in reshaped long_subset_data:", nrow(long_subset_data), "\n")

## Rows in reshaped long_subset_data: 53234

# Create matrices for betas and ses
betas_mat <- dcast(long_subset_data, SNP ~ trait, value.var = "beta", fill = NA)
ses_mat <- dcast(long_subset_data, SNP ~ trait, value.var = "se", fill = NA)
snps <- betas_mat$SNP

# Convert matrices to the appropriate format for hyprcoloc
rownames(betas_mat) <- betas_mat$SNP
rownames(ses_mat) <- ses_mat$SNP

betas_mat <- as.matrix(betas_mat[, -1, with = FALSE])
ses_mat <- as.matrix(ses_mat[, -1, with = FALSE])

# Restore the SNP column as rownames
rownames(betas_mat) <- snps
rownames(ses_mat) <- snps

# Confirm rownames are correctly set
head(rownames(betas_mat))

## [1] "rs1000364" "rs1002778" "rs1003323" "rs1003767" "rs1003768" "rs1003769"

cat("Dimensions of beta matrix:", dim(betas_mat), "\n")

## Dimensions of beta matrix: 4850 11

cat("Dimensions of SE matrix:", dim(ses_mat), "\n")

## Dimensions of SE matrix: 4850 11

# Check for rows with missing values
rows_with_na_betas <- apply(betas_mat, 1, function(row) any(is.na(row)))
rows_with_na_ses <- apply(ses_mat, 1, function(row) any(is.na(row)))

# Filter matrices to keep only rows with no missing values
betas_mat <- betas_mat[!rows_with_na_betas, , drop = FALSE]
ses_mat <- ses_mat[!rows_with_na_ses, , drop = FALSE]

cat("Dimensions of beta matrix:", dim(betas_mat), "\n")

## Dimensions of beta matrix: 4747 11

cat("Dimensions of SE matrix:", dim(ses_mat), "\n")

## Dimensions of SE matrix: 4747 11

# Run hyprcoloc
cat("Running hyprcoloc...\n")

## Running hyprcoloc...

betas <- betas_mat
ses <- ses_mat
traits_names <- colnames(betas_mat)
snp_ids <- rownames(betas_mat)

# Uniform priors = FALSE
results_unipriorsFALSE <- hyprcoloc(betas, ses, trait.names = traits_names, snp.id = snp_ids, snpscores = TRUE,uniform.priors=FALSE)

# Static table with kable
kable(
  results_unipriorsFALSE$results,
  caption = "HyprColoc Results with Uniform Priors = FALSE",
  align = 'c'
)

HyprColoc Results with Uniform Priors = FALSE
iteration	traits	posterior_prob	regional_prob	candidate_snp	posterior_explained_by_snp	dropped_trait
1	1-(1-enyl-palmitoyl)-2-palmitoyl-GPC (P-16:0/16:0), PCSK9, cholesterol, glycosyl ceramide (d18:1/20:0, d16:1/22:0), glycosyl-N-palmitoyl-sphingosine (d18:1/16:0), hydroxypalmitoyl sphingomyelin (d18:1/16:0(OH)), palmitoyl dihydrosphingomyelin (d18:0/16:0), palmitoyl sphingomyelin (d18:1/16:0), palmitoyl-sphingosine-phosphoethanolamine (d18:1/16:0), sphingomyelin (d18:1/24:1, d18:2/24:0)	0.9994	1	rs11591147	1	NA

# Uniform priors = TRUE
results_unipriorsTRUE <- hyprcoloc(betas, ses, trait.names = traits_names, snp.id = snp_ids, snpscores = TRUE, uniform.priors=TRUE)

kable(
  results_unipriorsTRUE$results,
  caption = "HyprColoc Results with Uniform Priors = TRUE",
  align = 'c'
)

HyprColoc Results with Uniform Priors = TRUE
iteration	traits	posterior_prob	regional_prob	candidate_snp	posterior_explained_by_snp	dropped_trait
1	1-(1-enyl-palmitoyl)-2-palmitoyl-GPC (P-16:0/16:0), PCSK9, cholesterol, glycosyl ceramide (d18:1/20:0, d16:1/22:0), glycosyl-N-palmitoyl-sphingosine (d18:1/16:0), hydroxypalmitoyl sphingomyelin (d18:1/16:0(OH)), palmitoyl dihydrosphingomyelin (d18:0/16:0), palmitoyl sphingomyelin (d18:1/16:0), palmitoyl-sphingosine-phosphoethanolamine (d18:1/16:0), sphingomyelin (d18:1/24:1, d18:2/24:0)	1	1	rs11591147	1	NA

# Without assuming uniform priors and with reg.thresh and align.thresh
results_uniform_priors_FALSE_reg.thresh.95_align.thresh.95 <- hyprcoloc(betas, ses, trait.names = traits_names, snp.id = snp_ids, uniform.priors = FALSE,snpscores = TRUE,
  reg.thresh = 0.95, align.thresh = 0.95
)

kable(
  results_uniform_priors_FALSE_reg.thresh.95_align.thresh.95$results,
  caption = "HyprColoc Results with Uniform Priors = FALSE and regional and alignment thresholds = 0.95",
  align = 'c'
)

HyprColoc Results with Uniform Priors = FALSE and regional and alignment thresholds = 0.95
iteration	traits	posterior_prob	regional_prob	candidate_snp	posterior_explained_by_snp	dropped_trait
1	1-(1-enyl-palmitoyl)-2-palmitoyl-GPC (P-16:0/16:0), PCSK9, cholesterol, glycosyl ceramide (d18:1/20:0, d16:1/22:0), glycosyl-N-palmitoyl-sphingosine (d18:1/16:0), hydroxypalmitoyl sphingomyelin (d18:1/16:0(OH)), palmitoyl dihydrosphingomyelin (d18:0/16:0), palmitoyl sphingomyelin (d18:1/16:0), palmitoyl-sphingosine-phosphoethanolamine (d18:1/16:0), sphingomyelin (d18:1/24:1, d18:2/24:0)	0.9994	1	rs11591147	1	NA

# Save the results as tabs in Excel
# Define file path for the output Excel workbook
output_file <- "/Users/charleenadams/pcsk9_hyprcoloc_corrected_results_top10_smallest_summary_with_snpscores.xlsx"

# Initialize a new workbook
wb <- createWorkbook()

# Function to add a result set to the workbook
add_results_to_workbook <- function(workbook, result_data, sheet_name) {
  # Extract the results data frame
  results_df <- result_data$results
  
  # Create a worksheet and write the data
  addWorksheet(workbook, sheet_name)
  writeData(workbook, sheet_name, results_df)
  
  # Apply styling
  header_style <- createStyle(fontSize = 12, fontColour = "#FFFFFF", halign = "CENTER", 
                              fgFill = "#4F81BD", textDecoration = "Bold")
  addStyle(workbook, sheet_name, header_style, rows = 1, cols = 1:ncol(results_df), gridExpand = TRUE)
}

# Function to add SNP scores to the workbook with rsIDs
add_snpscores_to_workbook <- function(workbook, snp_scores, rsIDs, sheet_name) {
  # Convert SNP scores to a data table and include rsIDs
  snp_scores_df <- as.data.table(snp_scores)
  colnames(snp_scores_df) <- paste0("SNP_Score_Set_", seq_len(ncol(snp_scores_df)))
  snp_scores_df$rsID <- rsIDs
  
  # Reorder columns to put rsID first
  setcolorder(snp_scores_df, c("rsID", setdiff(names(snp_scores_df), "rsID")))
  
  # Create a worksheet and write the data
  addWorksheet(workbook, sheet_name)
  writeData(workbook, sheet_name, snp_scores_df)
  
  # Apply styling
  header_style <- createStyle(fontSize = 12, fontColour = "#FFFFFF", halign = "CENTER", 
                              fgFill = "#4F81BD", textDecoration = "Bold")
  addStyle(workbook, sheet_name, header_style, rows = 1, cols = 1:ncol(snp_scores_df), gridExpand = TRUE)
}

# Add results and SNP scores for each analysis
# Ensure rsIDs are available as rownames in the matrices
rsIDs <- rownames(betas_mat)

# results_unipriorsFALSE
add_results_to_workbook(wb, results_unipriorsFALSE, "Results_Unipriors_FALSE")
add_snpscores_to_workbook(wb, results_unipriorsFALSE$snpscores, rsIDs, "SNP_Scores_Unipriors_FALSE")

# results_unipriorsTRUE
add_results_to_workbook(wb, results_unipriorsTRUE, "Results_Unipriors_TRUE")
add_snpscores_to_workbook(wb, results_unipriorsTRUE$snpscores, rsIDs, "SNP_Scores_Unipriors_TRUE")

# results_uniform_priors_FALSE_reg.thresh.95_align.thresh.95
add_results_to_workbook(wb, results_uniform_priors_FALSE_reg.thresh.95_align.thresh.95, 
                        "Unipriors_FALSE_Thresh_95")
add_snpscores_to_workbook(wb, results_uniform_priors_FALSE_reg.thresh.95_align.thresh.95$snpscores, 
                          rsIDs, "SNP_Scores_Thresh_95")

# Save the workbook
saveWorkbook(wb, output_file, overwrite = TRUE)

# Confirm completion
cat("Results and SNP scores saved to:", output_file, "\n")

## Results and SNP scores saved to: /Users/charleenadams/pcsk9_hyprcoloc_corrected_results_top10_smallest_summary_with_snpscores.xlsx

1.7 Summary tables: easier on the eyes

Here I provide James Staley’s (Hyprcoloc package author) example results. Their paper for the package had examined CHD. Their top putatively causal SNP (rs11591147) is the same as ours:

1.7.1 Staley’s table

Iteration	Traits	Posterior Probability	Regional Probability	Candidate SNP	Posterior Explained by SNP	Dropped Trait
1	T1, T2, T3, T4, T5	1.0000	1	rs11591147	1.0000	NA
2	T6, T7, T8	0.9164	1	rs12117612	0.4197	NA
3	T9, T10	0.9018	1	rs7524677	0.0763	NA

1.7.2 Our PCSK9 and 184 metabolites (uniform priors = FALSE)

Iteration	Traits	Posterior Prob	Regional Prob	Candidate SNP	Posterior Explained by SNP
1	1,2-dipalmitoyl-GPC (16:0/16:0), 1-(1-enyl-palmitoyl)-2-oleoyl-GPC (P-16:0/18:1), 1-(1-enyl-palmitoyl)-2-palmitoyl-GPC (P-16:0/16:0), 1-palmitoyl-2-stearoyl-GPC (16:0/18:0), 1-palmityl-2-palmitoyl-GPC (O-16:0/16:0), 1-stearoyl-2-oleoyl-GPI (18:0/18:1), 1-stearoyl-GPC (18:0), N-palmitoyl-sphingosine (d18:1/16:0), PCSK9, X - 23641, arachidoylcarnitine (C20), behenoyl sphingomyelin (d18:1/22:0), cerotoylcarnitine (C26), cholesterol, glycosyl ceramide (d18:1/20:0, d16:1/22:0), glycosyl ceramide (d18:2/24:1, d18:1/24:2), glycosyl-N-palmitoyl-sphingosine (d18:1/16:0), hydroxypalmitoyl sphingomyelin (d18:1/16:0(OH)), lactosyl-N-nervonoyl-sphingosine (d18:1/24:1), lignoceroyl sphingomyelin (d18:1/24:0), myristoyl dihydrosphingomyelin (d18:0/14:0), palmitoyl dihydrosphingomyelin (d18:0/16:0), palmitoyl sphingomyelin (d18:1/16:0), palmitoyl-sphingosine-phosphoethanolamine (d18:1/16:0), sphingomyelin (d17:1/16:0, d18:1/15:0, d16:1/17:0), sphingomyelin (d18:0/20:0, d16:0/22:0), sphingomyelin (d18:1/20:0, d16:1/22:0), sphingomyelin (d18:1/24:1, d18:2/24:0), sphingomyelin (d18:1/25:0, d19:0/24:1, d20:1/23:0, d19:1/24:0), sphingomyelin (d18:2/16:0, d18:1/16:1)	0.4593	0.6101	rs11591147	1

1.7.3 Our PCSK9 and 184 metabolites (uniform priors = TRUE)

Iteration	Traits	Posterior Prob	Regional Prob	Candidate SNP	Posterior Explained by SNP
1	1,2-dipalmitoyl-GPC (16:0/16:0), 1-(1-enyl-palmitoyl)-2-oleoyl-GPC (P-16:0/18:1), 1-(1-enyl-palmitoyl)-2-palmitoyl-GPC (P-16:0/16:0), 1-palmitoyl-2-stearoyl-GPC (16:0/18:0), 1-palmityl-2-palmitoyl-GPC (O-16:0/16:0), 1-stearoyl-2-oleoyl-GPI (18:0/18:1), 1-stearoyl-GPC (18:0), N-palmitoyl-sphingosine (d18:1/16:0), PCSK9, X - 23641, arachidoylcarnitine (C20), behenoyl sphingomyelin (d18:1/22:0), cerotoylcarnitine (C26), cholesterol, glycosyl ceramide (d18:1/20:0, d16:1/22:0), glycosyl ceramide (d18:2/24:1, d18:1/24:2), glycosyl-N-palmitoyl-sphingosine (d18:1/16:0), hydroxypalmitoyl sphingomyelin (d18:1/16:0(OH)), lactosyl-N-nervonoyl-sphingosine (d18:1/24:1), lignoceroyl sphingomyelin (d18:1/24:0), myristoyl dihydrosphingomyelin (d18:0/14:0), palmitoyl dihydrosphingomyelin (d18:0/16:0), palmitoyl sphingomyelin (d18:1/16:0), palmitoyl-sphingosine-phosphoethanolamine (d18:1/16:0), sphingomyelin (d17:1/16:0, d18:1/15:0, d16:1/17:0), sphingomyelin (d18:0/20:0, d16:0/22:0), sphingomyelin (d18:1/20:0, d16:1/22:0), sphingomyelin (d18:1/24:1, d18:2/24:0), sphingomyelin (d18:1/25:0, d19:0/24:1, d20:1/23:0, d19:1/24:0), sphingomyelin (d18:2/16:0, d18:1/16:1)	0.7712	0.7792	rs11591147	1
2	(S)-3-hydroxybutyrylcarnitine, 3-hydroxydecanoate, 3-hydroxydecanoylcarnitine, 3-hydroxyhexanoate, 3-hydroxylaurate, 3-hydroxyoctanoate, X - 15469, X - 16397, X - 18921, X - 21353, X - 22519, cis-4-decenoate (10:1n6), dodecadienoate (12:2), tetradecadienoate (14:2)	0.7638	0.7864	rs1616691	1
3	2,3-dihydroxy-5-methylthio-4-pentenoate (DMTPA), 2-palmitoyl-GPE (16:0), X - 15503, adipoylcarnitine (C6-DC), dimethylarginine (SDMA + ADMA), kynurenine, pseudouridine	0.7604	0.7791	rs141072298	0.8198
4	N-acetylkynurenine (2), cis-4-decenoylcarnitine (C10:1), decanoylcarnitine (C10), hexanoylcarnitine (C6), octanoylcarnitine (C8)	0.8019	0.8254	rs148238018	0.9829

1.7.4 Our PCSK9 and 184 metabolites (uniform priors = FALSE; regional and alignment thresholds = 0.95)

Iteration	Traits	Posterior Prob	Regional Prob	Candidate SNP	Posterior Explained by SNP
1	1,2-dipalmitoyl-GPC (16:0/16:0), 1-(1-enyl-palmitoyl)-2-oleoyl-GPC (P-16:0/18:1), 1-(1-enyl-palmitoyl)-2-palmitoyl-GPC (P-16:0/16:0), 1-palmitoyl-2-stearoyl-GPC (16:0/18:0), 1-palmityl-2-palmitoyl-GPC (O-16:0/16:0), N-palmitoyl-sphingosine (d18:1/16:0), PCSK9, behenoyl sphingomyelin (d18:1/22:0), cholesterol, glycosyl ceramide (d18:1/20:0, d16:1/22:0), glycosyl ceramide (d18:2/24:1, d18:1/24:2), glycosyl-N-palmitoyl-sphingosine (d18:1/16:0), hydroxypalmitoyl sphingomyelin (d18:1/16:0(OH)), lactosyl-N-nervonoyl-sphingosine (d18:1/24:1), lignoceroyl sphingomyelin (d18:1/24:0), myristoyl dihydrosphingomyelin (d18:0/14:0), palmitoyl dihydrosphingomyelin (d18:0/16:0), palmitoyl sphingomyelin (d18:1/16:0), palmitoyl-sphingosine-phosphoethanolamine (d18:1/16:0), sphingomyelin (d17:1/16:0, d18:1/15:0, d16:1/17:0), sphingomyelin (d18:0/20:0, d16:0/22:0), sphingomyelin (d18:1/20:0, d16:1/22:0), sphingomyelin (d18:1/24:1, d18:2/24:0), sphingomyelin (d18:2/16:0, d18:1/16:1)	0.9687	0.9808	rs11591147	1
2	(S)-3-hydroxybutyrylcarnitine, 3-hydroxydecanoate, 3-hydroxydecanoylcarnitine, 3-hydroxyhexanoate, 3-hydroxylaurate, 3-hydroxyoctanoate, X - 15469, X - 21353, X - 22519	0.9832	0.9884	rs1616691	1

1.7.5 Our PCSK9 and 10 metabolites (uniform priors = FALSE)

Iteration	Traits	Posterior Prob	Regional Prob	Candidate SNP	Posterior Explained by SNP
1	1-(1-enyl-palmitoyl)-2-palmitoyl-GPC (P-16:0/16:0), PCSK9, cholesterol, glycosyl ceramide (d18:1/20:0, d16:1/22:0), glycosyl-N-palmitoyl-sphingosine (d18:1/16:0), hydroxypalmitoyl sphingomyelin (d18:1/16:0(OH)), palmitoyl dihydrosphingomyelin (d18:0/16:0), palmitoyl sphingomyelin (d18:1/16:0), palmitoyl-sphingosine-phosphoethanolamine (d18:1/16:0), sphingomyelin (d18:1/24:1, d18:2/24:0)	0.9994	1	rs11591147	1

1.7.6 Our PCSK9 and 10 metabolites (uniform priors = TRUE)

Iteration	Traits	Posterior Prob	Regional Prob	Candidate SNP	Posterior Explained by SNP
1	1-(1-enyl-palmitoyl)-2-palmitoyl-GPC (P-16:0/16:0), PCSK9, cholesterol, glycosyl ceramide (d18:1/20:0, d16:1/22:0), glycosyl-N-palmitoyl-sphingosine (d18:1/16:0), hydroxypalmitoyl sphingomyelin (d18:1/16:0(OH)), palmitoyl dihydrosphingomyelin (d18:0/16:0), palmitoyl sphingomyelin (d18:1/16:0), palmitoyl-sphingosine-phosphoethanolamine (d18:1/16:0), sphingomyelin (d18:1/24:1, d18:2/24:0)	1	1	rs11591147	1

1.7.7 Our PCSK9 and 10 metabolites (uniform priors = FALSE; regional and alignment thresholds = 0.95)

Iteration	Traits	Posterior Prob	Regional Prob	Candidate SNP	Posterior Explained by SNP
1	1-(1-enyl-palmitoyl)-2-palmitoyl-GPC (P-16:0/16:0), PCSK9, cholesterol, glycosyl ceramide (d18:1/20:0, d16:1/22:0), glycosyl-N-palmitoyl-sphingosine (d18:1/16:0), hydroxypalmitoyl sphingomyelin (d18:1/16:0(OH)), palmitoyl dihydrosphingomyelin (d18:0/16:0), palmitoyl sphingomyelin (d18:1/16:0), palmitoyl-sphingosine-phosphoethanolamine (d18:1/16:0), sphingomyelin (d18:1/24:1, d18:2/24:0)	0.9994	1	rs11591147	1

Colocalization Across Clinical and Molecular Phenotypes

שׁוֹשַׁנָּה🌹

January 15, 2025