#!/bin/bash 

# Produced by Signature Science, LLC 4/7/2022

# Reproducibly run hap-ibd and IBIS with parameter sweeps to reproduce data and figures
# Requirements: 
# IBIS in path
# hap-ibd.jar copied to local directory
# all simulations unzipped in local directory (e.g. tar -xvfpz phased_simulations.tar.gz)


input_map="plink.allchr.GRCh37.map"
maps=genetic_map_hg19.txt
threads=10
for sim in *.vcf.gz
do
 prefix=${sim%%.vcf.gz}
 # Plink output required for ibis
 plink --vcf $sim --make-bed --out $prefix.plink  \
  && add-map-plink.pl $prefix.plink.bim $maps >$prefix.new.bim \
  && mv $prefix.new.bim $prefix.plink.bim
 
 # IBIS
 for ibis_er in 0.01  0.1 0.2
 do
  for ibis_mt in 200 100 50 10 5
  do
   for ibis_ml in 2 5 7
   do
    ibis -ibd2 -mt2 $ibis_mt -mL2 $ibis_ml -er2 $ibis_er -o ibis.$prefix.er-$ibis_er-$ibis_mt-$ibis_ml -printCoef -t $threads -b $prefix.plink
   done
  done
 done
# Hap-ibd
 for hap_gap in 1000 5000 10000
 do 
  for hap_mark in 100 50 10 2
  do
   for hap_len in {7,2,1}
   do
   hap_prefix=hap_gap-$hap_gap.hap_mark-$hap_mark.hap_len-$hap_len
   java -jar hap-ibd.jar \
        gt=$sim \
        map=$input_map \
        out=hap-ibd.${prefix}.${hap_prefix} \
        nthreads=$threads \
        max-gap=$hap_gap \
        min-output=$hap_len \
        min-markers=$hap_mark
   done
  done
 done

done
 


