Skip to content
Snippets Groups Projects
Commit 274f44a3 authored by Brandi Cantarel's avatar Brandi Cantarel
Browse files

Merge branch 'master' of git.biohpc.swmed.edu:BICF/Astrocyte/process_scripts

parents 6d29a48f 72c9c152
Branches
Tags
No related merge requests found
...@@ -67,7 +67,7 @@ bwa mem -M -t $NPROC -R "@RG\tID:${read_group}\tLB:tx\tPL:illumina\tPU:barcode\t ...@@ -67,7 +67,7 @@ bwa mem -M -t $NPROC -R "@RG\tID:${read_group}\tLB:tx\tPL:illumina\tPU:barcode\t
if [[ $umi == 'umi' ]] && [[ -f "${index_path}/genome.fa.alt" ]] if [[ $umi == 'umi' ]] && [[ -f "${index_path}/genome.fa.alt" ]]
then then
k8 ${testexe}/bwa-postalt.js -p tmphla ${index_path}/genome.fa.alt out.sam | python ${baseDir}/add_umi_sam.py -s - -o output.unsort.bam k8 ${testexe}/bwa-postalt.js -p tmphla ${index_path}/genome.fa.alt out.sam | python ${baseDir}/add_umi_sam.py -s - -o output.unsort.bam
elif [[ "${index_path}/genome.fa.alt" ]] elif [[ -f "${index_path}/genome.fa.alt" ]]
then then
k8 ${testexe}/bwa-postalt.js -p tmphla ${index_path}/genome.fa.alt out.sam| samtools view -1 - > output.unsort.bam k8 ${testexe}/bwa-postalt.js -p tmphla ${index_path}/genome.fa.alt out.sam| samtools view -1 - > output.unsort.bam
elif [[ $umi == 'umi' ]] elif [[ $umi == 'umi' ]]
......
...@@ -10,12 +10,13 @@ usage() { ...@@ -10,12 +10,13 @@ usage() {
exit 1 exit 1
} }
OPTIND=1 # Reset OPTIND OPTIND=1 # Reset OPTIND
while getopts :r:p:v:h opt while getopts :r:p:v:sh opt
do do
case $opt in case $opt in
p) pair_id=$OPTARG;; p) pair_id=$OPTARG;;
v) vcf=$OPTARG;; v) vcf=$OPTARG;;
r) index_path=$OPTARG;; r) index_path=$OPTARG;;
s) skipnorm=1;;
h) usage;; h) usage;;
esac esac
done done
...@@ -40,4 +41,9 @@ perl $baseDir\/uniform_vcf_gt.pl $pair_id $vcf ...@@ -40,4 +41,9 @@ perl $baseDir\/uniform_vcf_gt.pl $pair_id $vcf
bgzip -f ${pair_id}.uniform.vcf bgzip -f ${pair_id}.uniform.vcf
j=${pair_id}.uniform.vcf.gz j=${pair_id}.uniform.vcf.gz
tabix -f $j tabix -f $j
bcftools norm --fasta-ref $reffa -m - -Oz $j -o ${pair_id}.norm.vcf.gz if [[ skipnorm==1 ]]
then
cp $j ${pair_id}.norm.vcf.gz
else
bcftools norm --fasta-ref $reffa -m - -Oz $j -o ${pair_id}.norm.vcf.gz
fi
...@@ -79,7 +79,7 @@ while (my $line = <VCF>) { ...@@ -79,7 +79,7 @@ while (my $line = <VCF>) {
if ($hash{SVTYPE} eq 'DUP:TANDEM') { if ($hash{SVTYPE} eq 'DUP:TANDEM') {
print DUP join("\t",$chrom,$pos,$id,$ref,$alt,$score,$filter,$annot,$newformat,@newgts),"\n"; print DUP join("\t",$chrom,$pos,$id,$ref,$alt,$score,$filter,$annot,$newformat,@newgts),"\n";
}elsif ($hash{SVTYPE} eq 'DEL' || $hash{SVTYPE} eq 'INS') { }elsif ($hash{SVTYPE} eq 'DEL' || $hash{SVTYPE} eq 'INS') {
if (abs($hash{SVLEN}) < 50) { if (abs($hash{SVLEN}) < 20) {
print SI join("\t",$chrom,$pos,$id,$ref,$alt,$score,$filter,$annot,$newformat,@newgts),"\n"; print SI join("\t",$chrom,$pos,$id,$ref,$alt,$score,$filter,$annot,$newformat,@newgts),"\n";
}else { }else {
$newalt = "<".$hash{SVTYPE}.">"; $newalt = "<".$hash{SVTYPE}.">";
......
...@@ -104,7 +104,10 @@ then ...@@ -104,7 +104,10 @@ then
else else
/project/shared/bicf_workflow_ref/seqprg/svaba/bin/svaba run -p $NPROC -G ${reffa} -t ${sbam} -a ${pair_id} /project/shared/bicf_workflow_ref/seqprg/svaba/bin/svaba run -p $NPROC -G ${reffa} -t ${sbam} -a ${pair_id}
fi fi
java -Xmx10g -jar $SNPEFF_HOME/snpEff.jar -no-intergenic -lof -c $SNPEFF_HOME/snpEff.config ${snpeffgeno} ${pair_id}.svaba.unfiltered.somatic.sv.vcf | bgzip > ${pair_id}.svaba.vcf.gz vcf-concat ${pair_id}.svaba.unfiltered*sv.vcf | vcf-sort -t temp > svaba.unfiltered.vcf
bash $baseDir/norm_annot.sh -r ${index_path} -p svaba -v svaba.unfiltered.vcf -s
java -Xmx10g -jar $SNPEFF_HOME/snpEff.jar -no-intergenic -lof -c $SNPEFF_HOME/snpEff.config ${snpeffgeno} svaba.norm.vcf | bgzip > svaba.vcf.gz
fi fi
if [[ $method == 'lumpy' ]] if [[ $method == 'lumpy' ]]
......
...@@ -42,13 +42,17 @@ while (my $line = <VCF>) { ...@@ -42,13 +42,17 @@ while (my $line = <VCF>) {
foreach my $i (0..$#deschead) { foreach my $i (0..$#deschead) {
$gtdata{$deschead[$i]} = $gtinfo[$i]; $gtdata{$deschead[$i]} = $gtinfo[$i];
} }
if ($gtdata{AD}){ if ($gtdata{AD} =~ m/\d+,\d+/){
($gtdata{RO},@alts) = split(/,/,$gtdata{AD}); ($gtdata{RO},@alts) = split(/,/,$gtdata{AD});
$gtdata{AO} = join(",",@alts); $gtdata{AO} = join(",",@alts);
$gtdata{DP} = $gtdata{RO}; $gtdata{DP} = $gtdata{RO};
foreach (@alts) { foreach (@alts) {
$gtdata{DP} += $_; $gtdata{DP} += $_;
} }
} elsif ($gtdata{AD} =~ m/^\d+$/ && $gtdata{DP}){
$gtdata{AO} = $gtdata{AD};
$gtdata{RO} = $gtdata{DP} - $gtdata{AD};
$gtdata{AD} = join(',',$gtdata{RO},$gtdata{AO});
} elsif (exists $gtdata{NR} && exists $gtdata{NV}) { } elsif (exists $gtdata{NR} && exists $gtdata{NV}) {
$gtdata{DP} = $gtdata{NR}; $gtdata{DP} = $gtdata{NR};
$gtdata{AO} = $gtdata{NV}; $gtdata{AO} = $gtdata{NV};
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment