diff --git a/diff_exp/geneabundance.sh b/diff_exp/geneabundance.sh index 1baa572f82cf920d242ada960e3a7d86f1c2274f..e763ef7c10a6e3ec2f595fd1deb97f99d435b1c7 100644 --- a/diff_exp/geneabundance.sh +++ b/diff_exp/geneabundance.sh @@ -42,5 +42,5 @@ module load subread/1.6.1 stringtie/1.3.2d-intel featureCounts -s $stranded -M --fraction -J --ignoreDup -T $SLURM_CPUS_ON_NODE -p -g gene_name -a ${gtf} -o ${pair_id}.cts ${sbam} mkdir ${pair_id}_stringtie cd ${pair_id}_stringtie -stringtie ../${sbam} -p $SLURM_CPUS_ON_NODE -G ../${gtf} -B -e -o denovo.gtf -A ../${pair_id}.fpkm.txt +stringtie ../${sbam} -p $SLURM_CPUS_ON_NODE -G ${gtf} -B -e -o denovo.gtf -A ../${pair_id}.fpkm.txt diff --git a/genect_rnaseq/concat_ctsum.pl b/genect_rnaseq/concat_ctsum.pl new file mode 100644 index 0000000000000000000000000000000000000000..b6d83ef58fb5f9324d81ddc4cbc33dc6bc214fd1 --- /dev/null +++ b/genect_rnaseq/concat_ctsum.pl @@ -0,0 +1,45 @@ +#!/usr/bin/perl -w +#merge_tables.pl + +use Getopt::Long qw(:config no_ignore_case no_auto_abbrev); +use File::Basename; + +my %opt = (); +my $results = GetOptions (\%opt,'outdir|o=s','help|h'); + +$opt{outdir} = './' unless ($opt{outdir}); + +my @files = @ARGV; +foreach $file (@files) { + chomp($file); + open IN, "<$file" or die $!; + $fname = basename($file); + my $sample = (split(/\./,$fname))[0]; + $samps{$sample} = 1; + my $head = <IN>; + chomp($head); + my @colnames = split(/\t/,$head); + while (my $line = <IN>) { + chomp($line); + my @row = split(/\t/,$line); + my $gene = $row[0]; + my $ct = $row[-1]; + $cts{$gene}{$sample} = $ct; + } + close IN; +} + +open OUT, ">$opt{outdir}\/countFeature.txt" or die $!; +my @samples = sort {$a cmp $b} keys %samps; +print OUT join("\t",'FeatureCtType',@samples),"\n"; +foreach $gene (keys %cts) { + my @line; + foreach $sample (@samples) { + unless ($cts{$gene}{$sample}) { + $cts{$gene}{$sample} = 0; + } + push @line, $cts{$gene}{$sample}; + } + print OUT join("\t",$gene,@line),"\n"; +} +close OUT; diff --git a/genect_rnaseq/statanal.sh b/genect_rnaseq/statanal.sh new file mode 100644 index 0000000000000000000000000000000000000000..e221f1b8bd7b0fa72964cbf07a515065cef81988 --- /dev/null +++ b/genect_rnaseq/statanal.sh @@ -0,0 +1,43 @@ +#!/bin/bash +#rnaseqalign.sh + +usage() { + echo "-h Help documentation for hisat.sh" + echo "-d --DEA" + echo "Example: bash statanal.sh" + exit 1 +} +OPTIND=1 # Reset OPTIND +while getopts :d:h opt +do + case $opt in + d) dea=$OPTARG;; + h) usage;; + esac +done + +shift $(($OPTIND -1)) + +# Check for mandatory options +if [[ -z $SLURM_CPUS_ON_NODE ]] +then + SLURM_CPUS_ON_NODE=1 +fi +source /etc/profile.d/modules.sh + +perl $baseDir/scripts/concat_cts.pl -o ./ *.cts +perl $baseDir/scripts/concat_fpkm.pl -o ./ *.fpkm.txt +perl $baseDir/scripts/concat_ctsum.pl -o ./ *.cts.summary +cp design.txt design.shiny.txt +cp geneset.gmt geneset.shiny.gmt + +if [[ $dea == 'skip' ]] +then + touch empty.png + touch bg.rda +else + module load R/3.2.1-intel + Rscript $baseDir/scripts/dea.R + Rscript $baseDir/scripts/build_ballgown.R *_stringtie + perl $baseDir/scripts/concat_edgeR.pl *.edgeR.txt +fi