<jupyter_start><jupyter_text>Perform variant quality controlIn this notebook, we perform variant quality control using [PLINK2](https://www.cog-genomics.org/plink/2.0/).Note that this work is part of a larger project to [Demonstrate the Potential for Pooled Analysis of All of Us and UK Biobank Genomic Data](https://github.com/all-of-us/ukb-cross-analysis-demo-project). Specifically this is for the portion of the project that is the **pooled** analysis.  Setup     Cloud Environment: This notebook was written for use on the All of Us Workbench.            Use "Recommended Environment" General Analysis which creates compute type Standard VM with reasonable defaults for CPU, RAM and disk.        This notebook takes about 2 hours to run. Recommend that it is run in the background via run_notebook_in_the_background.<jupyter_code>from datetime import datetime
import os
import time<jupyter_output><empty_output><jupyter_text>Setup plink2https://www.cog-genomics.org/plink/2.0/<jupyter_code>%%bash

##### plink 2 install
PLINK_VERSION=2.3.Alpha
PLINK_ZIP_PATH=/tmp/plink-$PLINK_VERSION.zip
curl -L -o $PLINK_ZIP_PATH https://s3.amazonaws.com/plink2-assets/alpha2/plink2_linux_x86_64.zip
mkdir -p /tmp/plink2/
unzip -o $PLINK_ZIP_PATH -d /tmp/plink2/
!/tmp/plink2/plink2 --version # --help<jupyter_output><empty_output><jupyter_text>Define constants<jupyter_code># Papermill parameters. See https://papermill.readthedocs.io/en/latest/usage-parameterize.html

#---[ Inputs ]---
# The BGEN file created via aou_workbench_pooled_analyses/05_write_pooled_bgen.ipynb.
REMOTE_MERGED_BGEN = 'gs://fc-secure-e53e4a44-7fe2-42b7-89b7-01aae1e399f7/data/pooled/geno/20220215/aou-alpha3-ukb-chr1-chr22.bgen'
# The sample file created via aou_workbench_pooled_analyses/05_write_pooled_bgen.ipynb.
REMOTE_MERGED_BGEN_SAMPLE = 'gs://fc-secure-e53e4a44-7fe2-42b7-89b7-01aae1e399f7/data/pooled/geno/20220215/aou-alpha3-ukb-chr1-chr22.sample'
# This TSV was created via notebook aou_workbench_pooled_analyses/02_merge_lipids_phenotypes.ipynb
REMOTE_PHENOTYPE_SAMPLE_IDS = 'gs://fc-secure-e53e4a44-7fe2-42b7-89b7-01aae1e399f7/data/pooled/pheno/20220310/aou_alpha3_ukb_lipids_ids.tsv'

#---[ Outputs ]---
# Create a timestamp for a folder of results generated today.
DATESTAMP = time.strftime('%Y%m%d')
OUTPUT_FILENAME_PREFIX = 'aou_alpha3_ukb_lipids'
OUTPUT_FOLDER = f'{os.getenv("WORKSPACE_BUCKET")}/data/pooled/variant-qc/{DATESTAMP}/'
LOCAL_MERGED_BGEN = os.path.basename(REMOTE_MERGED_BGEN)
LOCAL_MERGED_BGEN_SAMPLE = os.path.basename(REMOTE_MERGED_BGEN_SAMPLE)
LOCAL_PHENOTYPE_SAMPLE_IDS = os.path.basename(REMOTE_PHENOTYPE_SAMPLE_IDS)<jupyter_output><empty_output><jupyter_text>Copy data locally<jupyter_code>!gsutil cp -n {REMOTE_MERGED_BGEN} {REMOTE_MERGED_BGEN_SAMPLE} .
!gsutil cp {REMOTE_PHENOTYPE_SAMPLE_IDS} .<jupyter_output><empty_output><jupyter_text>Variant QC for regenie step 1 via plink2Use plink2 to perform [variant QC appropriate for step 1 of regenie](https://rgcgithub.github.io/regenie/recommendations/exclusion-files), obtaining a subset of high quality variants roughly equal to the number of samples.<jupyter_code>!/tmp/plink2/plink2 \
  --bgen {LOCAL_MERGED_BGEN} ref-first \
  --sample {LOCAL_MERGED_BGEN_SAMPLE} \
  --chr 1-22 \
  --keep {LOCAL_PHENOTYPE_SAMPLE_IDS} \
  --geno 0.1 \
  --mind 0.1 \
  --mac 100 \
  --hwe 1e-15 \
  --write-snplist \
  --write-samples \
  --no-id-header \
  --out {OUTPUT_FILENAME_PREFIX}_step1QC_plink
!ls -lth . | head
!head {OUTPUT_FILENAME_PREFIX}_step1QC_plink.id
!wc -l {OUTPUT_FILENAME_PREFIX}_step1QC_plink.id
!head {OUTPUT_FILENAME_PREFIX}_step1QC_plink.snplist
!wc -l {OUTPUT_FILENAME_PREFIX}_step1QC_plink.snplist
!ls {OUTPUT_FILENAME_PREFIX}_step1QC_plink*
!gsutil -m cp {OUTPUT_FILENAME_PREFIX}_step1QC_plink* {OUTPUT_FOLDER}<jupyter_output><empty_output><jupyter_text>Variant QC for regenie step 2 via plink2Use plink2 to perform [variant QC appropriate for step 2 of regenie](https://rgcgithub.github.io/regenie/recommendations/exclusion-files), obtaining a subset of high quality variants.<jupyter_code>!/tmp/plink2/plink2 \
  --bgen {LOCAL_MERGED_BGEN} ref-first \
  --sample {LOCAL_MERGED_BGEN_SAMPLE} \
  --chr 1-22 \
  --keep {LOCAL_PHENOTYPE_SAMPLE_IDS} \
  --geno 0.1 \
  --mind 0.1 \
  --hwe 1e-15 \
  --write-snplist \
  --write-samples \
  --no-id-header \
  --out {OUTPUT_FILENAME_PREFIX}_step2QC_plink

# This differs from step 1 QC in that the following parameter was removed
#    --mac 100 \
!ls -lth {OUTPUT_FILENAME_PREFIX}*
!wc -l {OUTPUT_FILENAME_PREFIX}_step2QC_plink.id
!wc -l {OUTPUT_FILENAME_PREFIX}_step2QC_plink.snplist
!ls {OUTPUT_FILENAME_PREFIX}_step2QC_plink*
!gsutil -m cp {OUTPUT_FILENAME_PREFIX}_step2QC_plink* {OUTPUT_FOLDER}
!gsutil ls {OUTPUT_FOLDER}<jupyter_output><empty_output><jupyter_text>Provenance<jupyter_code>%%bash

date
%%bash

pip3 freeze<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.10974358974358975, 'fraction_numerical': 0.03569230769230769, 'mean_word_length': 5.283505154639175, 'pattern_counts': {'":': 0, '<': 29, '<?xml version=': 0, '>': 29, 'https://': 7, 'lorem ipsum': 0, 'www.': 2, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '7', 'dir': 'jupyter-structured-clean-dedup', 'id': '354722', 'license': 'BSD-3-Clause', 'n_tokens_mistral': 2057, 'n_tokens_neox': 1878, 'n_words': 354, 'path': 'aou_workbench_pooled_analyses/06_pooled_variant_qc.ipynb', 'repo_name': 'all-of-us/ukb-cross-analysis-demo-project', 'size': '4875'}