-
Notifications
You must be signed in to change notification settings - Fork 7
New issue
Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.
By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.
Already on GitHub? Sign in to your account
Unify host genome generation #182
base: main
Are you sure you want to change the base?
Conversation
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Great stuff, thanks @morsecodist!
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Looks good just a couple of nits
runtime { | ||
docker: docker_image_id | ||
cpu: cpu | ||
memory: "240G" |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Do we want this to be either GiB or GB?
String docker_image_id | ||
} | ||
|
||
call ensure_gz as genome_fasta { |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Are we getting the results of this function anywhere?
TMPDIR=${TMPDIR:-/tmp} | ||
|
||
if [ "~{nucleotide_type}" == "dna" ]; then | ||
>&2 minimap2 -x map-ont -d '~{genome_name}_~{nucleotide_type}.mmi' "~{fasta}" |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Not sure why, but nucleotide_type
is not getting filled in with dna
or rna
but genome_name
is.
/hisat2/hisat2_extract_splice_sites.py <(pigz -dc '~{transcripts_gtf_gz}') > "$TMPDIR/genome.ss" & pid=$! | ||
/hisat2/hisat2_extract_exons.py <(pigz -dc '~{transcripts_gtf_gz}') > "$TMPDIR/genome.exon" | ||
wait $pid | ||
>&2 /hisat2/hisat2-build -p 16 \ |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
@mlin @morsecodist I'm having trouble with this build when the only transcripts_gtf_gz
is the ERCC.gtf. The hisat2-build command fails with Error: Encountered exception: 'Nongraph exception'
The ERCC.gtf is needed for STAR to properly run.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Can we say if the file TMPDIR/genome.ss
don't bother running with --exon "$TMPDIR/genome.exon" --ss "$TMPDIR/genome.ss"
?
* fix ercc gtf in index generation * pigz to cat
set -euxo pipefail | ||
TMPDIR=${TMPDIR:-/tmp} | ||
|
||
mkdir -p "$TMPDIR"'/bt2/~{genome_name}' |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
There's a discrepancy here with the old host filtering where the name of the folder of the old host filtering is
"~{genome_name}.bowtie2"
This is causing inconsistency between the two and making it so that the new bowtie2 indices are not back compatible. We could create a symlink here before tarring to make it back-compatible. Otherwise we'd have to regenerate either the modern or old host filtering indices
* fastp * fastp single * bowtie2 run * hisat2 run * dedup run * run subsample * run kallisto * adjust index tar filenames * polishing * polishing * count reads in each step * Create host_filter_indexing.wdl * boost fastp complexity threshold * output fastp report * build fastp from our fork with SDUST complexity filtering * use fastp --sdust_complexity_filter * bump * bump * tune * stub the remaining step descriptions * wire to tests * and auto_benchmark * fixup tests * fixup tests * fixup tests * fixup tests * fixup tests * fixup tests * add back in picard CollectInsertSizeMetrics * picard step description * host_filter_2022.wdl => host_filter.wdl * polish * restore fastqs_0 and fastqs_1 to minimize collateral changes * add minimap2 index build * picard_insert_metrics.txt * amr/run.wdl workaround * index multiple transcripts_fasta_gz * make gtf optional * allow uncompressed genome fasta * allow uncompressed genome fasta * allow uncompressed genome fasta * bump minimap2 memory * bump minimap2 memory * step descriptions -- first draft * add indexing driver & draft readme * include invocations in step descriptions * rebase amr fix * load card_json * run kallisto every time * fix amr wdl * fix short-read-mngs rebase weirdness * add final things * [modernized host filter] add ERCC and gene-level outputs to kallisto (#175) The kallisto step gains two new derivative output files: * `ERCC_counts.tsv`: Estimated read counts for the ERCC sequences only (two-column TSV: ERCC_id, est_counts) * `gene_abundance.tsv`: gene-level est_counts and tpm, computed by summing over all transcripts for each gene * (and `abundance.tsv` is renamed to `transcript_abundance.tsv`) To get the `gene_abundance.tsv` we need a new input `gtf_gz`, the Ensembl GTF file for the host species that will tell it how to map the transcript IDs in `transcript_abundance.tsv` onto gene IDs for the roll-up. The input is optional and if absent then the `gene_abundance.tsv` output is omitted too. Note: docker image update needed to install & upgrade some dependencies. * load card_json explicitly * add ~ * fix host_filter unit tests * fix host_filter unit tests * bowtie2: sort by read name for better reproducibility * update minimap2 indexing invocation * add chelonia_mydas, drosophila_melanogaster, gray_whale, pea-aphid * copy-paste {bowtie2,hisat2}_human_filter to support pipeline viz * allow kallisto nonzero exit * rename modern host filtering inputs/outputs and create a 1-1 mapping between inputs/outputs * fix lint issue * rename reads_in_count to input_read_count * auto_benchmark updates * fix test_RunCZIDDedup_safe_csv * rename kallisto output files * update mosquitos with several Culicidae * add files to wdl output for pipeline viz compatibility * convert headers in descriptions to bolded text * delete host_filter_indexing since it's subsumed in #182 * fix glob patterns in read counting * Revert "fix glob patterns in read counting" This reverts commit aeb234f. * [Bug] fix count expansion for single file short-read-mngs (#216) * fix bowtie2 counts for single file * fix extra expansions * relieve hisat2 dependency * single sample hisat2 * fix hisat2 * fix dockerfile for hisat2 --------- Co-authored-by: Omar Valenzuela <51972068+ovalenzuela19@users.noreply.github.com> * Remove AMR changes that are a WIP from modern host filtering branch (#219) * Revert "output gene id in primary output file (#209)" This reverts commit 2d9ff56. * Revert "Output non host reads and non host contigs for AMR (#205)" This reverts commit 9de3fc2. * tune hisat2 memory usage (#223) * Legacy Host Filter initial commit (#224) * legacy-host-filter-inital-commit * linting * add stage io map * remove stage io map swp file * Revert "Remove AMR changes that are a WIP from modern host filtering branch (#219)" (#226) This reverts commit 227a489. --------- Co-authored-by: Mike Lin <mlin@Mikes-MacBook-Pro.local> Co-authored-by: Omar Valenzuela <ovalenzuela@chanzuckerberg.com> Co-authored-by: Omar Valenzuela <51972068+ovalenzuela19@users.noreply.github.com> Co-authored-by: rzlim08 <37033997+rzlim08@users.noreply.github.com>
No description provided.