metabat won't execute successfully

[ycsong]

Hello,

I am able to submit the jobs to the node succesfully, but I am having some issues with metabat.

First a quick question. For metabat, I see that there needs to be two types of input:

• A gzipped fasta file named contigs.fasta.gz placed in the assemblies/<sample_code> folder.
• A 'cov' placed in metabat/<sample_code> folder

What is this 'cov' file, or is it a folder? I am assuming that it is a coverage stats table from assemblies step that looks like:

#ID	Avg_fold	Length	Ref_GC	Covered_percent	Covered_bases	Plus_reads	Minus_reads	Read_GC	Median_fold	Std_Dev
CFS1-BTM1_c1	20.2561	2487	0.5690	99.3567	2471	166	170	0.5730	21	7.14
CFS1-BTM3_c1	22.2687	2021	0.6081	100.0000	2021	150	151	0.6106	24	9.39
CFS1-BTM4_c1	8.8488	2017	0.5786	100.0000	2017	58	60	0.5788	10	3.81
CFS1-BTM2_c1	31.9714	1995	0.5278	98.7970	1971	211	215	0.5032	42	21.06

It would be good if this could be confirmed.

My second question relates to setting up the conda environment. I saw it in the log that the following is trying to be executed:

set +u;source activate envs/metagem;set -u

I tried executing this code outside of metagem and then I got the following error message:

EnvironmentLocationNotFound: Not a conda environment: /home/bioinfo_tools/bin/metaGEM/workflow/envs/metagem

What would be a possible source of the issue?

Thank you in advance,

Young

[ycsong]

Update: I looked at this issue further and realized that I need to edit the "env" section of the config.yaml, such that whatever specified there is actually a conda environment. This appears to have solved the problem.

I do have a question: Is it possible to specify which data that I want to work with? For example, if I have subfolders A and B in the dataset folder, could I customize files or commands such that my run targets the B subfolder?

[franciscozorrilla]

Update: I looked at this issue further and realized that I need to edit the "env" section of the config.yaml, such that whatever specified there is actually a conda environment. This appears to have solved the problem.

Glad you figured that one out!

What is this 'cov' file, or is it a folder? I am assuming that it is a coverage stats table from assemblies step that looks like:

As I have recommended in our previous discussion, please have a look at the Snakefile rules to understand what is going on:

metaGEM/workflow/Snakefile

Lines 390 to 492 in 57cfd60

	rule crossMapSeries:
	input:
	contigs = rules.megahit.output,
	reads = f'{config["path"]["root"]}/{config["folder"]["qfiltered"]}'
	output:
	concoct = directory(f'{config["path"]["root"]}/{config["folder"]["concoct"]}/{{IDs}}/cov'),
	metabat = directory(f'{config["path"]["root"]}/{config["folder"]["metabat"]}/{{IDs}}/cov'),
	maxbin = directory(f'{config["path"]["root"]}/{config["folder"]["maxbin"]}/{{IDs}}/cov')
	benchmark:
	f'{config["path"]["root"]}/{config["folder"]["benchmarks"]}/{{IDs}}.crossMapSeries.benchmark.txt'
	message:
	"""
	Cross map in seies:
	Use this approach to provide all 3 binning tools with cross-sample coverage information.
	Will likely provide superior binning results, but may no be feasible for datasets with
	many large samples such as the tara oceans dataset.
	"""
	shell:
	"""
	# Activate metagem environment
	set +u;source activate {config[envs][metagem]};set -u;
	# Create output folders
	mkdir -p {output.concoct}
	mkdir -p {output.metabat}
	mkdir -p {output.maxbin}
	# Make job specific scratch dir
	idvar=$(echo $(basename $(dirname {output.concoct})))
	echo -e "\nCreating temporary directory {config[path][scratch]}/{config[folder][crossMap]}/${{idvar}} ... "
	mkdir -p {config[path][scratch]}/{config[folder][crossMap]}/${{idvar}}
	# Move into scratch dir
	cd {config[path][scratch]}/{config[folder][crossMap]}/${{idvar}}
	# Copy files
	cp {input.contigs} .
	# Define the focal sample ID, fsample:
	# The one sample's assembly that all other samples' read will be mapped against in a for loop
	fsampleID=$(echo $(basename $(dirname {input.contigs})))
	echo -e "\nFocal sample: $fsampleID ... "
	echo "Renaming and unzipping assembly ... "
	mv $(basename {input.contigs}) $(echo $fsampleID|sed 's/$/.fa.gz/g')
	gunzip $(echo $fsampleID|sed 's/$/.fa.gz/g')
	echo -e "\nIndexing assembly ... "
	bwa index $fsampleID.fa
	for folder in {input.reads}/*;do
	id=$(basename $folder)
	echo -e "\nCopying sample $id to be mapped against the focal sample $fsampleID ..."
	cp $folder/*.gz .
	# Maybe I should be piping the lines below to reduce I/O ?
	echo -e "\nMapping sample to assembly ... "
	bwa mem -t {config[cores][crossMap]} $fsampleID.fa *.fastq.gz > $id.sam
	echo -e "\nConverting SAM to BAM with samtools view ... "
	samtools view -@ {config[cores][crossMap]} -Sb $id.sam > $id.bam
	echo -e "\nSorting BAM file with samtools sort ... "
	samtools sort -@ {config[cores][crossMap]} -o $id.sort $id.bam
	echo -e "\nRunning jgi_summarize_bam_contig_depths script to generate contig abundance/depth file for maxbin2 input ... "
	jgi_summarize_bam_contig_depths --outputDepth $id.depth $id.sort
	echo -e "\nMoving depth file to sample $fsampleID maxbin2 folder ... "
	mv $id.depth {output.maxbin}
	echo -e "\nIndexing sorted BAM file with samtools index for CONCOCT input table generation ... "
	samtools index $id.sort
	echo -e "\nRemoving temporary files ... "
	rm *.fastq.gz *.sam *.bam
	done
	nSamples=$(ls {input.reads}|wc -l)
	echo -e "\nDone mapping focal sample $fsampleID agains $nSamples samples in dataset folder."
	echo -e "\nRunning jgi_summarize_bam_contig_depths for all sorted bam files to generate metabat2 input ... "
	jgi_summarize_bam_contig_depths --outputDepth $id.all.depth *.sort
	echo -e "\nMoving input file $id.all.depth to $fsampleID metabat2 folder... "
	mv $id.all.depth {output.metabat}
	echo -e "Done. \nCutting up contigs to 10kbp chunks (default), not to be used for mapping!"
	cut_up_fasta.py -c {config[params][cutfasta]} -o 0 -m $fsampleID.fa -b assembly_c10k.bed > assembly_c10k.fa
	echo -e "\nSummarizing sorted and indexed BAM files with concoct_coverage_table.py to generate CONCOCT input table ... "
	concoct_coverage_table.py assembly_c10k.bed *.sort > coverage_table.tsv
	echo -e "\nMoving CONCOCT input table to $fsampleID concoct folder"
	mv coverage_table.tsv {output.concoct}
	echo -e "\nRemoving intermediate sorted bam files ... "
	rm *.sort
	"""

In summary, this rule will generate the mapping information needed for each individual binner, and store it under the corresponding binner/sample_ID/cov/ subfolder. Then, each binning job you submit will look into the appropriate folder in order to fetch the coverage information needed for binning :)

I suggest you run the crossMapSeries rule above to see what the output files and folders look like. You may have to delete any empty /cov directories for the jobs to run, in case Snakemake didnt do that automatically when your last jobs failed.

I do have a question: Is it possible to specify which data that I want to work with? For example, if I have subfolders A and B in the dataset folder, could I customize files or commands such that my run targets the B subfolder?

This is more of a general Snakemake question, so I recommend reading the documentation for this. As a quick workaround, you could simply temporarily relocate the /dataset subfolders you are not interested out of that folder, then the Snakefile wont be able to find those samples to submit jobs.

[ycsong]

Thanks for your help Francisco. I think at this point, I can experiment with few things in Snakefile and other config documents to get things running (I've already did a couple of test runs in the server, and it seems to run well). If I do have more q's in the future, I'll post them in a separate discussion.

Thanks again,

Young