Skip to content

Commit

Permalink
add read count script
Browse files Browse the repository at this point in the history
  • Loading branch information
sanjaynagi committed Aug 12, 2024
1 parent dadc5fc commit e9332f0
Showing 1 changed file with 38 additions and 0 deletions.
38 changes: 38 additions & 0 deletions workflow/scripts/countreads.sh
Original file line number Diff line number Diff line change
@@ -0,0 +1,38 @@
#!/bin/bash

# Set reference files
GFF_FILE="resources/reference/Anopheles-gambiae-PEST_BASEFEATURES_AgamP4.12.gff3"

# Sample names
SAMPLES=("Kisumu1" "Kisumu2" "Kisumu3" "Kisumu4" "G28-BusiaSurvivors1" "G28-BusiaSurvivors2" "G28-BusiaSurvivors3" "G28-BusiaSurvivors4")

# Loop through each sample and run the tools
for SAMPLE in "${SAMPLES[@]}"
do
BAM_FILE="results/alignments/${SAMPLE}.hisat2.bam"
R1_FASTQ="resources/reads/${SAMPLE}_1.fastq"
R2_FASTQ="resources/reads/${SAMPLE}_2.fastq"

FEATURECOUNTS_OUTPUT="${SAMPLE}_featureCounts_output.txt"
HTSEQCOUNT_OUTPUT="${SAMPLE}_htseq_counts_output.txt"
KALLISTO_OUTPUT="${SAMPLE}_kallisto_output"

# Run featureCounts
echo "Running featureCounts for ${SAMPLE}..."
featureCounts -a $GFF_FILE -o $FEATURECOUNTS_OUTPUT -p -B -t exon -g Parent -C $BAM_FILE
if [ $? -ne 0 ]; then
echo "featureCounts failed for ${SAMPLE}"
exit 1
fi

# Run htseq-count
echo "Running htseq-count for ${SAMPLE}..."
htseq-count -f bam -r pos -s no -t exon -i Parent $BAM_FILE $GFF_FILE > $HTSEQCOUNT_OUTPUT
if [ $? -ne 0 ]; then
echo "htseq-count failed for ${SAMPLE}"
exit 1
fi

done

echo "All tools have completed successfully for all samples."

0 comments on commit e9332f0

Please sign in to comment.