nievergeltlab’s gists

nievergeltlab / make_locuszoom.qsub

Created May 1, 2017 18:54

Filter METAL outputs for LocusZoom - Splits results by chromosome

	#!/bin/bash
	#PBS -V

	while getopts l:d: option
	do
	case "${option}"
	in
	l) snplocations=${OPTARG};;
	d) lzresults_list=${OPTARG};;
	esac

nievergeltlab / merge_discordance.r

Created April 19, 2017 20:55

Descriptive statistics for genotype merge discordances

	#Merge data in PLINK, using merge mode 6 or 7
	./plink --bfile YEHUDA --bmerge YEHUDA.bed YEHUDA.bim YEHUDA.fam --merge-mode 7 --out yehude-merge

	R
	setwd('F:/rutgers_2')
	dat <- read.table('yehude-merge.diff', header=T,nr=800000,stringsAsFactors=F)
	library(plyr)
	#Determine general amount of disagreement for each SNP. Ones that have especially high disagreement may be badly genotyped
	quantile(table(dat$SNP))
	#Plot it

nievergeltlab / convert_dosage_to_gmmat_format.sh

Created April 10, 2017 16:00

Logistic Mixed Models with GMMAT

	#Script to convert dosage data to GMMAT format. Requires also a list of usuable SNPs as input, but this list does not need to be filtered.

	#Create a valid SNPlist
	info=0.9
	maf=0.01
	zcat /home/maihofer/vets/qc/imputation/distribution/vets_eur_analysis_run2/daner_vets_eur_analysis_run2.gz \| awk -v info=$info -v maf=$maf '{if (($8 > info) && ($6 > maf) && ($7 > maf) && ($6 < 1-maf) && ($7 < 1-maf) && ($11 != "NA")) print $2}' > european_valid_snps.snplist
	echo "SNP" > snp.txt

	cat snp.txt european_valid_snps.snplist \| LC_ALL=C sort -k 1b,1 > european_valid_snps.sorted.snplist

nievergeltlab / tar_wildcards.sh

Created April 8, 2017 23:17

Extract files from tar with wildcards

tar xvzf dnhs_qc_v2_mar7_2017.tgz --wildcards '*-qc.fam' --wildcards '*.header'

nievergeltlab / beta_logscale_function.r

Last active April 13, 2018 17:10

Implementation of formula converting beta from a linear regression onto the log odds scale EFFICIENT COMPUTATION WITH A LINEAR MIXED MODEL ON LARGE-SCALE DATA SETSWITH APPLICATIONS TO GENETIC STUDIES Author(s): Matti Pirinen, Peter Donnelly and Chri

	R

	new_beta_se <- function(x ,phi=0.153)
	{
	#Input should be a vector consisting of beta, p, maf. phi should be a value for the proportion of cases
	B=x[1]
	P=x[2]
	theta=x[3]

	B2=B/(phi(1 - phi) + .5(1-2phi)(1-2theta)B - (0.084 + 0.9phi(1 - 2phi)theta(1-theta))/(phi(1-phi))*B^2)

nievergeltlab / extract_snp_from_dosages

Created March 24, 2017 18:23

Algorithm: Find which file the SNP is in, then zgrep the SNP out. This script will output a .dosage file for each snp in the list. to combine them into a file, it's of course just a matter of using "cat"

	#Requires a list of SNPs written into snplist.txt

	dosage_directory=write_path_to_dosage_files_here
	for snp in $(cat snplist.txt)
	do
	filename=$(grep -w -m1 -l $snp "$dosage_directory"/*.map \| sed 's/.map/.gz/g' )
	zgrep -w -m1 $snp $filename > "$snp".dosage
	done

nievergeltlab / analyze_dosages.pbs

Last active April 4, 2017 19:59

Analyze dosages data from Ricopili

	#PBS -lnodes=1
	#PBS -lwalltime=0:05:00

	#!/bin/bash

	while getopts a:b:o:p:d: option
	do
	case "${option}"
	in
	a) phenotype=${OPTARG};;

nievergeltlab / ibd_compare

Created March 22, 2017 22:08

2 file IBD, if strand mismatches and allele mismatches

	bimfile=pts_meg2_mix_am-qc.bim
	grep -P "A\tT" $bimfile > ambiguous_snps.txt
	grep -P "T\tA" $bimfile >> ambiguous_snps.txt
	grep -P "C\tG" $bimfile >> ambiguous_snps.txt
	grep -P "G\tC" $bimfile >> ambiguous_snps.txt


	#filer out ambiguous
	$plink_location --bfile pts_meg2_mix_am-qc --exclude ambiguous_snps.txt --make-bed --out pts_meg2_mix_am-qc-ambig

nievergeltlab / name_chromosome_chunk.sh

Created March 13, 2017 16:18

Get chromosome number and chunk number from file name. For ricopili.

	chunknum=$(echo $chunk \| sed 's/.$chr[0-9]_[0-9]_[0-9]$.*/\1/')
	chrnum=$(echo $chunk \| sed 's/.$chr[0-9]$.*/\1/')

nievergeltlab / manhattan_plot.r

Last active May 1, 2017 18:51

Manhattan Plots V2 of the plotter will color SNPs within a set kb window of any hits

	args <- commandArgs(trailingOnly = TRUE)
	scriptloc <- args[1]
	results <- args[2]
	outfile <- args[3]
	snpcol <- args[4]
	chrcol <- args[5]
	bpcol <- args[6]
	pcol <- args[7]

	chunknum=$(echo $chunk \| sed 's/.\(chr[0-9]_[0-9]_[0-9]\).*/\1/')
	chrnum=$(echo $chunk \| sed 's/.\(chr[0-9]\).*/\1/')

Adam X. Maihofer nievergeltlab