From 1b3892e5395d5f5aa2da2e2c4e6ca8580ed40a91 Mon Sep 17 00:00:00 2001
From: Brandi Cantarel <brandi.cantarel@utsouthwestern.edu>
Date: Mon, 13 Jul 2020 17:22:16 -0500
Subject: [PATCH] qc for research
---
alignment/bamqc.sh | 16 +++++++++++++---
alignment/sequenceqc_dna.pl | 4 +++-
2 files changed, 16 insertions(+), 4 deletions(-)
diff --git a/alignment/bamqc.sh b/alignment/bamqc.sh
index 79db6ad..80a32a2 100644
--- a/alignment/bamqc.sh
+++ b/alignment/bamqc.sh
@@ -35,6 +35,11 @@ shift $(($OPTIND -1))
# usage
#fi
+if [[ -z $version ]]
+then
+ version='NA'
+fi
+
source /etc/profile.d/modules.sh
module load samtools/gcc/1.10 fastqc/0.11.8
samtools flagstat ${sbam} > ${pair_id}.flagstat.txt
@@ -55,7 +60,7 @@ fi
tmpdir=`pwd`
if [[ $nuctype == 'dna' ]]; then
module load bedtools/2.29.2 picard/2.10.3
- bedtools coverage -sorted -g ${index_path}/genomefile.txt -a ${bed} -b ${sbam} -hist > ${pair_id}.covhist.txt
+ bedtools coverage -a ${bed} -b ${sbam} -hist > ${pair_id}.covhist.txt
grep ^all ${pair_id}.covhist.txt > ${pair_id}.genomecov.txt
perl $baseDir/calculate_depthcov.pl ${pair_id}.covhist.txt
if [[ -z $skiplc ]]
@@ -63,13 +68,18 @@ if [[ $nuctype == 'dna' ]]; then
samtools view -@ $NPROC -b -L ${bed} -o ${pair_id}.ontarget.bam ${sbam}
samtools index -@ $NPROC ${pair_id}.ontarget.bam
samtools flagstat ${pair_id}.ontarget.bam > ${pair_id}.ontarget.flagstat.txt
- samtools view -@ $NPROC -b -q 1 ${sbam} | bedtools coverage -sorted -hist -g ${index_path}/genomefile.txt -b stdin -a ${bed} > ${pair_id}.mapqualcov.txt
+ samtools view -@ $NPROC -b -q 1 ${sbam} | bedtools coverage -hist -b stdin -a ${bed} > ${pair_id}.mapqualcov.txt
java -Xmx64g -Djava.io.tmpdir=${tmpdir} -XX:ParallelGCThreads=$NPROC -jar $PICARD/picard.jar EstimateLibraryComplexity BARCODE_TAG=RG I=${sbam} OUTPUT=${pair_id}.libcomplex.txt TMP_DIR=${tmpdir}
#java -Xmx64g -Djava.io.tmpdir=${tmpdir} -jar $PICARD/picard.jar CollectAlignmentSummaryMetrics R=${index_path}/genome.fa I=${pair_id}.ontarget.bam OUTPUT=${pair_id}.alignmentsummarymetrics.txt TMP_DIR=${tmpdir}
#samtools view -@ $NPROC ${sbam} | awk '{sum+=$5} END { print "Mean MAPQ =",sum/NR}' > ${pair_id}.meanmap.txt
fi
#java -Xmx64g -Djava.io.tmpdir=${tmpdir} -jar $PICARD/picard.jar CollectInsertSizeMetrics INPUT=${sbam} HISTOGRAM_FILE=${pair_id}.hist.ps REFERENCE_SEQUENCE=${index_path}/genome.fa OUTPUT=${pair_id}.hist.txt TMP_DIR=${tmpdir}
- perl $baseDir/sequenceqc_dna.pl -e ${version} -r $index_path ${pair_id}.genomecov.txt
+ if [[ $index_path/reference_info.pl ]]
+ then
+ perl $baseDir/sequenceqc_dna.pl -e ${version} -r $index_path ${pair_id}.genomecov.txt
+ else
+ touch ${pair_id}.genomecov.txt
+ fi
else
perl $baseDir/sequenceqc_rna.pl -e ${version} -r $index_path ${pair_id}.flagstat.txt
fi
diff --git a/alignment/sequenceqc_dna.pl b/alignment/sequenceqc_dna.pl
index 9d573b8..6634923 100755
--- a/alignment/sequenceqc_dna.pl
+++ b/alignment/sequenceqc_dna.pl
@@ -122,7 +122,9 @@ foreach $sfile (@statfiles) {
"Alignment_Status\t".$hash{status},"Alignment_Date\t".$hash{date},
"File_Owner\t".$hash{fileowner},"Workflow_Version\t".$gittag),"\n";
close OUT;
- system(qq{cat $opt{refdir}\/reference_info.txt >> $prefix\.sequence.stats.txt});
+ if (-e "$opt{refdir}\/reference_info.txt") {
+ system(qq{cat $opt{refdir}\/reference_info.txt >> $prefix\.sequence.stats.txt});
+ }
##### END separateFilesPerSample ######
}
--
GitLab