-
Notifications
You must be signed in to change notification settings - Fork 3
Expand file tree
/
Copy pathxContrastTwins.sh
More file actions
93 lines (74 loc) · 2.66 KB
/
xContrastTwins.sh
File metadata and controls
93 lines (74 loc) · 2.66 KB
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
#!/bin/bash
#$ -cwd -l mem=1G,time=:30: -N FilterTwin
#get arguments
while getopts v:t:n:p:f:a:b:u: opt; do
case "$opt" in
v) VcfFil="$OPTARG";;
t) TwnTab="$OPTARG";;
n) DirPre="$OPTARG";;
p) AddPrm="$OPTARG";;
f) FamNam="$OPTARG";;
a) Twin1="$OPTARG";;
b) Twin2="$OPTARG";;
u) Unfiltered="$OPTARG";;
#H) echo "$usage"; exit;;
esac
done
FiltScrDir="/home/local/ARCS/ads2202/scripts/Filtering_scripts"
VcfFil=`readlink -f $VcfFil`
if [[ $TwnTab ]]; then
FamNam=`cut -f 1 $TwnTab | head -n $SGE_TASK_ID | tail -n 1`
Twin1=`cut -f 2 $TwnTab | head -n $SGE_TASK_ID | tail -n 1`
Twin2=`cut -f 3 $TwnTab | head -n $SGE_TASK_ID | tail -n 1`
Unfiltered=`cut -f 4 $TwnTab | head -n $SGE_TASK_ID | tail -n 1`
fi
if [[ $FamNam == [0-9]* ]]; then FamNam=Fam$FamNam; fi
OutFil=$FamNam.ContrastTwins
DirNam=$FamNam
if [[ $DirPre ]]; then DirNam=$DirPre"_"$DirNam; fi
mkdir -p $DirNam
cd $DirNam
runCustomFilter(){
if [[ $Unfiltered ]]; then TwinComp=$TwinComp" --unfl $Unfiltered"; fi
if [[ $AddPrm ]]; then TwinComp=$TwinComp" $AddPrm"; fi
CMD="$FiltScrDir/ExmFilt.CustomGenotype.py -v $VcfFil $TwinComp"
echo $CMD
eval $CMD
}
TwinComp=" -o $OutFil.t1.temp --het $Twin1 --ref $Twin2"
runCustomFilter
TwinComp=" -o $OutFil.t2.temp --het $Twin1 --alt $Twin2"
runCustomFilter
TwinComp=" -o $OutFil.t3.temp --ref $Twin1 --het $Twin2"
runCustomFilter
TwinComp=" -o $OutFil.t4.temp --ref $Twin1 --alt $Twin2"
runCustomFilter
TwinComp=" -o $OutFil.t5.temp --alt $Twin1 --ref $Twin2"
runCustomFilter
TwinComp=" -o $OutFil.t6.temp --alt $Twin1 --het $Twin2"
runCustomFilter
R --vanilla <<RSCRIPT
options(stringsAsFactors=F)
dat <- read.delim("$OutFil.t1.temp.tsv")
for(i in 2:6){
temp <- read.delim(paste("$OutFil.t", i, ".temp.tsv", sep=""))
temp <- temp[,match(colnames(temp), colnames(dat))]
dat <- rbind(dat, temp)
}
write.table(dat, "$OutFil.tsv", row.names=F, quote=F, sep="\t")
RSCRIPT
LEN=`cat $OutFil.tsv | wc -l`
if [[ $LEN -gt 1 ]]; then
qsub $FiltScrDir/xAnnotateVariantTSV.sh -i $OutFil.tsv
fi
echo "Six contrasts performed on paired samples. Logs for each contrast:" > $OutFil.filtered.log
echo >> $OutFil.filtered.log
echo "-------------------------------------------------------------------------------------------------------------" >> $OutFil.filtered.log
echo >> $OutFil.filtered.log
for i in $OutFil.t*.temp.log; do
grep -vE "OutputName|Not" $i >> $OutFil.log
echo >> $OutFil.filtered.log
echo "-------------------------------------------------------------------------------------------------------------" >> $OutFil.filtered.log
echo >> $OutFil.filtered.log
done
rm -f $OutFil.*temp*