GEM script that makes individual scripts for each phenotype

Same script was also run with BMI in –covar-names for BMI-adjusted analysis

Make sure to run GEM v.1.4.3 on intel nodes with the build that has been optimized with fast C++ libraries. This will improve genome-wide runtime from 2+ days to ~1 hour.

Using robust standard errors to account for P-value inflation due to heteroskedasticity.


#Script that makes scripts for each phenotype to run indiv by chromosome.

phenotypes=("Alanine_aminotransferase_rint" "Albumin_rint" "Alkaline_phosphatase_rint" 
"Apolipoprotein_A_rint" "Apolipoprotein_B_rint" "Aspartate_aminotransferase_rint" 
"C_reactive_protein_rint" "Calcium_rint" "Cholesterol_rint" "Creatinine_rint" 
"Cystatin_C_rint" "Direct_bilirubin_rint" "Gamma_glutamyltransferase_rint"
"Glucose_rint" "HbA1c_rint" "HDL_cholesterol_rint" "IGF_1_rint" 
"LDL_direct_rint" "Lipoprotein_A_rint" "Oestradiol_rint" 
"Phosphate_rint" "Rheumatoid_factor_rint" "SHBG_rint" "Testosterone_rint" 
"Total_bilirubin_rint" "Total_protein_rint" "Triglycerides_rint" 
"Urate_rint" "Urea_rint" "Vitamin_D_rint" "eGFR_rint" 
"bioavailableTest_rint" "freeTest_rint")

genoindir=("/scratch/mf91122/BioxVeg/genotypeQC/GEM2")
outdir=("/scratch/mf91122/BioxVeg/gem5")
pheno=("/scratch/mf91122/BioxVeg/pheno/09112022/BioxVeg_pheno_09112022.csv")

for j in ${phenotypes[@]} 
    do

mkdir -p /work/kylab/mike/BioxVeg/GEM5/byPheno/$j


echo "#!/bin/bash
#SBATCH --partition=batch
#SBATCH --job-name=GEM-BioxVeg-"$j"
#SBATCH --nodes=1
#SBATCH --ntasks-per-node=8
#SBATCH --time=144:00:00
#SBATCH --mem=30000
#SBATCH --output=%x.%j.out
#SBATCH --error=%x.%j.err
#SBATCH --constraint=EDR
#SBATCH --array=1-22

i=\$SLURM_ARRAY_TASK_ID

cd /work/kylab/mike/BioxVeg/GEM5/byPheno/$j
ml GEM/1.4.3-intel-2020b


mkdir -p $outdir/$j

GEM \
--bgen $genoindir/chr\$i.bgen \
--sample $genoindir/chr\$i.sample \
--pheno-file $pheno \
--sampleid-name IID \
--pheno-name $j \
--covar-names zAge PCA1 PCA2 PCA3 PCA4 PCA5 \
PCA6 PCA7 PCA8 PCA9 PCA10 PCA11 PCA12 PCA13 \
PCA14 PCA15 PCA16 PCA17 PCA18 PCA19 PCA20 \
zTownsend zblood Sex Geno_batch \
AlcoholFreq PreviousSmoker CurrentSmoker \
--robust 1 \
--exposure-names Veg5yr1 \
--threads 8 \
--output-style meta \
--out $outdir/$j/chr\$i

" > /work/kylab/mike/BioxVeg/GEM5/byPheno/$j/GEM5-$j.sh

#Uncomment to run all those scripts
sbatch /work/kylab/mike/BioxVeg/GEM5/byPheno/$j/GEM5-$j.sh


done #end pheno loop