fix: statistics from GWAS

[mcanouil]

Fixes

• In R/run_eggla_gwas.R,
  • fix: update path to files in example.
  • fix: check for samples in VCF files.
  • fix: leading zeros in IDs in dataset were removed.
  • fix: use PLINK2 to compute HWE P-value, N, alleles frequencies, and missing rate.
  • fix: use PLINK2 to compute MACH_R2 imputation score.
  • fix: use vep_file as the top level parameter instead of vep.
  • fix: remove --mach-r2-filter PLINK2 option.
• In R/plot_aucs.R, R/plot_egg_aucs.R, and R/run_eggla_lmm.R,
  • fix: groupOnX argument of geom_quasirandom() is deprecated as of ggbeeswarm 0.7.1.
• In R/plot_slopes.R and R/plot_egg_slopes.R,
  • fix: Using size aesthetic for lines was deprecated.

Chores

• chore: split .devcontainer settings into a separate repositories, i.e., .vscode and .devcontainer.

Build

• build: upgrade PLINK2 and BCFtools to latest versions.
• build: update R packages to latest versions.

---

After all discussions, below is the current output from the GWAS.

[mcanouil]

I added --freq, --missing, and --hardy from PLINK2. This will now produce the following output:

Note that FDR, BONFERRONI, etc. columns are only there to allow individual cohorts to reuse "locally" the summary statistics. It's not necessarily for use in the meta-analysis.

@annihei Let me know, if the changes suit the group, see #97 for more details.

[mcanouil]

@annihei Is the final word to remove HWE test?

Does the other changes seem fine to you? If so, could you mark the different threads above using the "Resolve conversation" buttons?

[annihei]

@mcanouil Yes let's drop the HWE_P. I have been trying to double check the other changes, will comment on those.

I think we should also drop the --mach-r2-filter flag in all PLINK2 commands as this was not in analysis plan and we filter based on INFO in the meta-analysis team.

[mcanouil]

I think we should also drop the --mach-r2-filter flag in all PLINK2 commands as this was not in analysis plan and we filter based on INFO in the meta-analysis team.

Removed in aa18ffc323b45c6cb9c4fcc39ad681a5f82390aa

[mcanouil]

HWE commands/outputs removed in c7c09c69c431d33c8dc5304a94808431f8248a98

[mcanouil]

@annihei You should be able to install the version from this PR using the following command: remotes::install_github("mcanouil/eggla@fix/issue96").

Note that, the column names might not match exactly the analysis plan (I did not check that yet). Either, we update the analysis plan to match the R package (It would make more sense that way to me), or we rename the output columns within the R package.

[annihei]

OK, thanks. I'll try these.

[mcanouil]

For reference, the current GWAS results output

It seems OBS_CT and MISSING_CT does not count DS data but only hardcall genotypes

[mcanouil]

Now, the values checks out. Note that Call rate is based on hard call for now. Should it be computed on dosage? F_MISS_DOSAGE as described in https://www.cog-genomics.org/plink/2.0/formats#vmiss

[mcanouil]

To-do:

• [x] check allele frequencies it it matches with dosage.
• [x] rename columns as with previous output.
• [x] remove or rename CALL_RATE which represents m, in the current analysis output, the genotyped rate.
• [x] in analytical plan add https://www.cog-genomics.org/plink/2.0/assoc#glm_errcode for ERRCODE description.

[mcanouil]

[mcanouil]

After all discussions, below is the current output from the GWAS.