SNPTEST GWAS script

run_snptest.sh

#!/bin/bash

#PBS -N sleepgwas
#PBS -o sleepsnptest-output
#PBS -e sleepsnptest-error
#PBS -t 1-22
#PBS -l walltime=24:00:00
#PBS -l nodes=1:ppn=2
#PBS -S /bin/bash

set -e

echo "Running on ${HOSTNAME}"

if [ -n "${1}" ]; then
  echo "${1}"
  PBS_ARRAYID=${1}
fi

i=${PBS_ARRAYID}

module add apps/snptest.2.5.0


genfile="/panfs/panasas01/shared/alspac/studies/latest/alspac/genetic/variants/arrays/gwas/imputed/1000genomes/released/27Feb2015/data/genotypes/dosage/data_chr${i}"
samplefile="/panfs/panasas01/shared/alspac/studies/latest/alspac/genetic/variants/arrays/gwas/imputed/1000genomes/released/27Feb2015/data/data.sample"

# Change this file to point to where you want
outfile="${HOME}/results_chr${i}"

snptest_v2.5 \
	-data ${genfile} ${samplefile} \
	-missing_code -9 \
	-pheno plink_pheno \
	-cov_all \
	-use_raw_phenotypes \
	-frequentist 1 \
	-method em \
	-o ${outfile}

This is a PBS job submission script that will run a GWAS using SNPTEST on bluecrystal3, making 22 parallel jobs (one for each chromosome).

To use for your own GWAS just change the "samplefile" name and the "outfile" name (and the -N, -o and -e flags). Then in the frontend of bluecrystal3 simply run:

qsub run_snptest.sh

Alternatively, if you want to test that the script is running e.g. for chromosome 22, you can run

chmod 755 run_snptest.sh
./run_snptest.sh 22

and it will start to run the script with the chromosome variable set to 22. Avoid running in the front end for too long, i.e. just use this to test it's doing as you expect. To kill press ctrl+c.

Important thing is to make sure that the sample file is constructed correctly - use the .sample file in the imputation data directory as a template and this represents the correct order of individuals in the gen file.

Optional extra flags:

• -exclude_snps <snplist file>
• -exclude_samples <sample list file>

See https://mathgen.stats.ox.ac.uk/genetics_software/snptest/snptest.html for more details

When running batch jobs in bluecrystal3 you can see if the set of jobs is running like this:

qstat -u <username>

or you can see if specific jobs within the batch are running like this:

showq -u <username>