raw-lab / mercat2

MerCat2: python code for versatile k-mer counting and diversity estimation for database independent property analysis for metaome data

MerCat2: python code for versatile k-mer counter for database independent property analysis (DIPA) for omic analysis

Installing MerCat2

Option 1: Bioconda Installer

Install mamba using conda

conda activate base
conda install mamba

• NOTE: Make sure you install mamba in your base conda environment
• We have found that mamba is faster than conda for installing packages and creating environments. Using conda might fail to resolve dependencies.

Install MerCat2

• Available via Bioconda:

mamba create -n mercat2 -c conda-forge -c bioconda mercat2
conda activate mercat2

Option 2: Source Installer

• Clone mercat2 from github
• Run install_mercat2.sh to install all required dependencies
• This script creates a conda environment for you

git clone https://github.com/raw-lab/mercat2.git
cd mercat2
bash install_mercat2.sh
conda activate mercat2

Dependencies

MerCat2 runs on python version 3.9 and higher.

external dependencies

MerCat2 can run without external dependencies based on the options used.

Required dependencies:

• When a raw read .fastq file is given

  • Fastqc
  • Fastp

• For bacteria/archaea rich samples (-prod option)

  • Prodigal

• For eukaryote rich samples or general applications (-fgs option)

  • FragGeneScanRs

These are available through BioConda, except FragGeneScanRS, which is included in the MerCat2 distribution.

conda install -c bioconda fastqc fastp prodigal

Usage

usage: mercat2.py [-h] [-i I [I ...]] [-f F] -k K [-n N] [-c C] [-prod] [-fgs] [-s S] [-o O] [-replace] [-lowmem LOWMEM] [-skipclean] [-toupper] [-pca] [--version]

options:
  -h, --help      show this help message and exit
  -i I            path to input file(s)
  -f F            path to folder containing input files
  -k K            kmer length
  -n N            no of cores [auto detect]
  -c C            minimum kmer count [10]
  -prod           run Prodigal on fasta files
  -fgs            run FragGeneScanRS on fasta files
  -s S            Split into x MB files. [100]
  -o O            Output folder, default = 'mercat_results' in current directory
  -replace        Replace existing output directory [False]
  -lowmem LOWMEM  Flag to use incremental PCA when low memory is available. [auto]
  -skipclean      skip trimming of fastq files
  -toupper        convert all input sequences to uppercase
  -pca            create interactive PCA plot of the samples (minimum of 4 fasta files required)
  --version, -v   show the version number and exit

Mercat assumes the input file format based on the extension provided

• raw fastq file: ['.fastq', '.fq']

• nucleotide fasta: ['.fa', '.fna', '.ffn', '.fasta']

• amino acid fasta: ['.faa']

• It also accepts gzipped versions of these filetypes with the added '.gz' suffix

Usage examples

Run mercat2 on a protein file (protein fasta - '.faa')

mercat2.py -i file-name.faa -k 3 -c 10

Run mercat2 on a nucleotide file (nucleotide fasta - '.fa', '.fna', '.ffn', '.fasta')

mercat2.py -i file-name.fna -k 3 -n 8 -c 10

Run mercat2 on a nucleotide file raw data (nucleotide fastq - '.fastq')

mercat2.py -i file-name.fastq -k 3 -n 8 -c 10

Run on many samples within a folder

mercat2.py -f /path/to/input-folder -k 3 -n 8 -c 10

Run on sample with prodigal/FragGeneScanRS option (raw reads or nucleotide contigs - '.fa', '.fna', '.ffn', '.fasta', '.fastq')

mercat2.py -i /path/to/input-file -k 3 -n 8 -c 10 -prod

mercat2.py -i /path/to/input-file -k 3 -n 8 -c 10 -fgs

• the prodigal and FragGeneScanRS options run the k-mer counter on both contigs and produced amino acids

Outputs

• Results are stored in the output folder (default 'mercat_results' of the current working directory)
  • the 'report' folder contains an html report with interactive plotly figures
    • If at least 4 samples are provided a PCA plot will be included in the html report
  • the 'tsv' folder contains counts tables in tab separated format
    • if protein files are given, or the -prod option, a .tsv file is created for each sample containing k-mer count, pI, Molecular Weight, and Hydrophobicity metrics
    • if nucleotide files are given a .tsv file is created for each sample containing k-mer count and GC content
  • if .fastq raw reads files are used, a 'clean' folder is created with the clean fasta file.
  • if the -prod option is used, a 'prodigal' folder is created with the amino acid .faa and .gff files
  • if the -fgs option is used, a 'fgs' folder is created with the amino acid .faa file

Diversity estimation

Alpha and Beta diversity metrics provided by MerCat2 are experimental. We are currently working on the robustness of these measures.

Alpha diversity metrics provided:

• shannon
• simpson
• simpson_e
• goods_coverage
• fisher_alpha
• dominance
• chao1
• chao1_ci
• ace

Beta diversity metrics provided:

• euclidean
• cityblock
• braycurtis
• canberra
• chebyshev
• correlation
• cosine
• dice
• hamming
• jaccard
• mahalanobis
• manhattan (same as City Block in this case)
• matching
• minkowski
• rogerstanimoto
• russellrao
• seuclidean
• sokalmichener
• sokalsneath
• sqeuclidean
• yule

Notes on memory usage and speed

MerCat2 uses a substantial amount of memory when the k-mer is high.
Running MerCat2 on a personal computer using a k-mer length of ~4 should be OK. Total memory usage can be reduced using the Chunker feature (-s option), but keep in mind that in testing when the chunk size is too small (1MB) some of the least significant k-mers will get lost. This does not seem to affect the overall results, but it is something to keep in mind. Using the chunker and reducing the number of CPUs available (-n option) can help reduce memory requirements.

The speed of MerCat2 can be increased when more memory or computer nodes are available on a cluster and using a chunk size of about 100Mb.

Copyright

This is copyrighted by University of North Carolina at Charlotte, Jose L. Figueroa III, Andrew Redinbo, and Richard Allen White III. All rights reserved. DeGenPrime is a bioinformatic tool that can be distributed freely for academic use only. Please contact us for commerical use. The software is provided “as is” and the copyright owners or contributors are not liable for any direct, indirect, incidental, special, or consequential damages including but not limited to, procurement of goods or services, loss of use, data or profits arising in any way out of the use of this software.

Citing Mercat

If you are publishing results obtained using MerCat2, please cite:
Figueroa JL*, Redinbo A*, Panyala A, Colby S, Friesen M, Tiemann L, White III RA. 2024.
MerCat2: a versatile k-mer counter and diversity estimator for database-independent property analysis obtained from omics data
Bioinformatics Advances, vbae061 Bioinformatics Advances
*Co-first authors

BioRxiv pre-print
Figueroa JL, Panyala A, Colby S, Friesen M, Tiemann L, White III RA. 2022.
MerCat2: a versatile k-mer counter and diversity estimator for database-independent property analysis obtained from omics data.
bioRxiv

CONTACT

Please send all queries to Jose Luis Figueroa III
Dr. Richard Allen White III
Andrew Redinbo

Or open an issue

mercat2's Issues

Recommendations

Hi @raw-lab @raw937

Thanks for your previous answer regarding the difference between mercat and mercat2. I was wondering if you have any recommendations on:

1. Number of input reads to use (I have illumina reads and I've tried downsample 100k, 1M, 5M...etc)
2. kmer length (I've used k=21, which depending on input can impact speed of running)
3. I was thinking of using mercat2 to quickly calculate the shannon/simpsons index as a way to compare diversity between samples (running mercat individually on each sample). I was wondering if that can be used that way, as a way to compare between samples.

thanks again for your help!

Will

Difference between mercat and mercat2?

Hi

I was wondering what the main difference between mercat and mercat2 is. From the readme, it looks like the difference is the output, and the plots? Is that correct or are there any differences in the underlying methods.

thanks!

IndexError: list index out of range

Hi,

I would like to analyze the similarity of a set of genome sequence (Vibrio stains) with mercat2's approach.
I installed mercat2 with mamba as instructed in the repository (using Python 3.10.14), and runned as follows.
I noticed the error happen only when multi sample (2≥) were supplied to the program.
How should I deal with this ?

best,

=============================================================================
$ cp GCF*fna test/
$ mercat2.py -f test/ -k 3 -n 8 -c 10
.
.
...
GCF_024442235.1 .fna
GCF_028743655.1 .fna
GCF_025917725.1 .fna
GCF_003691505.1 .fna
GCF_009665315.1 .fna
GCF_003691545.1 .fna
GCF_001310575.2 .fna
GCF_022453845.1 .fna
GCF_002887655.1 .fna
GCF_002021755.1 .fna
Time to load 16 files: 8.51 seconds
Checking for large nucleotide files
Processing Nucleotides
Running Mercat2 using 18 cores
Time to count 3-mers: 1.71 seconds

Creating Nucleotide Graphs
/home/kazu/mambaforge/envs/mercat2/lib/python3.10/site-packages/sklearn/metrics/pairwise.py:2317: DataConversionWarning:

Data was converted to boolean for metric dice

/home/kazu/mambaforge/envs/mercat2/lib/python3.10/site-packages/sklearn/metrics/pairwise.py:2317: DataConversionWarning:

Data was converted to boolean for metric jaccard

Error with beta metric: Mahalanobis
Data must be symmetric and cannot contain NaNs.
/home/kazu/mambaforge/envs/mercat2/lib/python3.10/site-packages/sklearn/metrics/pairwise.py:2317: DataConversionWarning:

Data was converted to boolean for metric matching

/home/kazu/mambaforge/envs/mercat2/lib/python3.10/site-packages/sklearn/metrics/pairwise.py:2317: DataConversionWarning:

Data was converted to boolean for metric rogerstanimoto

/home/kazu/mambaforge/envs/mercat2/lib/python3.10/site-packages/sklearn/metrics/pairwise.py:2317: DataConversionWarning:

Data was converted to boolean for metric russellrao

/home/kazu/mambaforge/envs/mercat2/lib/python3.10/site-packages/sklearn/metrics/pairwise.py:2317: DataConversionWarning:

Data was converted to boolean for metric sokalmichener

/home/kazu/mambaforge/envs/mercat2/lib/python3.10/site-packages/sklearn/metrics/pairwise.py:2317: DataConversionWarning:

Data was converted to boolean for metric sokalsneath

/home/kazu/mambaforge/envs/mercat2/lib/python3.10/site-packages/sklearn/metrics/pairwise.py:2317: DataConversionWarning:

Data was converted to boolean for metric yule

Gathering Diversity Metrics
Traceback (most recent call last):
File "/home/kazu/mambaforge/envs/mercat2/bin/mercat2.py", line 508, in
mercat_main()
File "/home/kazu/mambaforge/envs/mercat2/bin/mercat2.py", line 499, in mercat_main
mercat2_report.merge_tsv(tomerge, outfile)
File "/home/kazu/mambaforge/envs/mercat2/lib/python3.10/site-packages/mercat2_lib/mercat2_report.py", line 128, in merge_tsv
kmer = sorted(kmers)[0]
IndexError: list index out of range

Lack of diversity report after the completion of the run

Hey,

I wanted to view some diversity metrics for my dataset but unfortunately, no file was generated containing the diversity metrics of the dataset. The program outputted a protein summary.txt file which contained numbers but with no description what the numbers meant and the html report wanted me download a .tsv file but I have not been able to successfully download the file. Any suggestion or help or clarification is welcome.

Thanks

ImportError: cannot import name 'ConstantInputWarning' from 'scipy.stats'

Error due to scikit-kbio:

ImportError: cannot import name 'ConstantInputWarning' from 'scipy.stats' (/home/smith/miniconda3/lib/python3.8/site-packages/scipy/stats/__init__.py)

This is likely related to recent installations (like numba) which might have modified versions of numpy and scipy, as previously the problem didn't occur.

Versions:
python : 3.8
scipy : 1.8.1
numpy : 1.22.4
skbio : 0.5.9

Related:
Error: cannot import name 'SpearmanRConstantInputWarning' from 'scipy.stats'
Cannot import name 'ConstantInputWarning' from 'scipy.stats' #4

Option for transposed table output

Thank you very much for this package, it's been a great help!

It would be very helpful to have an option to output transposed tables, with the kmers along the columns and the fasta files as rows. I'm currently using an awk script to transpose them, but when k and the number of genomes get larger, memory becomes an issue.

Cannot import name 'ConstantInputWarning' from 'scipy.stats'

Hello, this tool looks great! But when running mercat2 using demo data, I encountered the following error:

File "/storage/home/hcoda1/1/kbian8/.conda/envs/mercat2/lib/python3.10/site-packages/skbio/stats/distance/_mantel.py", line 16, in
from scipy.stats import ConstantInputWarning
ImportError: cannot import name 'ConstantInputWarning' from 'scipy.stats' (/storage/home/hcoda1/1/kbian8/.conda/envs/mercat2/lib/python3.10/site-packages/scipy/stats/init.py)

SciPy has been updated to 1.8.1 and all dependencies have been installed as well as their version matched.

No idea how to solve this. May I have your suggestions?

Failed to start dashboard / cannot import name 'packaging' from 'pkg_resources' / Unable to register worker with raylet. No such file or directory

I am trying to run mercat2 on SwissProt. In a little test invocation,

   mercat2.py -i uniprot_sprot.fasta -k1 -c1 -skipclean -o k1

I get the below errors.

I also tried starting a Ray cluster manually by

   ray start --head

but that gives me another error as below.

What am I doing wrong?

Many thanks,
David.

Error messages from mercat2 alone:

Starting MerCat2 v1.4.1 with k-mer 1 and 48 threads

2024-05-23 19:52:54,205 ERROR services.py:1207 -- Failed to start the dashboard , return code 1
2024-05-23 19:52:54,206 ERROR services.py:1232 -- Error should be written to 'dashboard.log' or 'dashboard.err'. We are printing the last 20 lines for you. See 'https://docs.ray.io/en/master/ray-observability/ray-logging.html#logging-directory-structure' to find where the log file is.
2024-05-23 19:52:54,206 ERROR services.py:1276 -- 
The last 20 lines of /tmp/ray/session_2024-05-23_19-52-52_585467_40696/logs/dashboard.log (it contains the error message from the dashboard): 
  File "/bi/home/dkreil/.conda/envs/mercat2/lib/python3.10/site-packages/ray/dashboard/dashboard.py", line 75, in run
    await self.dashboard_head.run()
  File "/bi/home/dkreil/.conda/envs/mercat2/lib/python3.10/site-packages/ray/dashboard/head.py", line 322, in run
    modules = self._load_modules(self._modules_to_load)
  File "/bi/home/dkreil/.conda/envs/mercat2/lib/python3.10/site-packages/ray/dashboard/head.py", line 219, in _load_modules
    head_cls_list = dashboard_utils.get_all_modules(DashboardHeadModule)
  File "/bi/home/dkreil/.conda/envs/mercat2/lib/python3.10/site-packages/ray/dashboard/utils.py", line 121, in get_all_modules
    importlib.import_module(name)
  File "/bi/home/dkreil/.conda/envs/mercat2/lib/python3.10/importlib/__init__.py", line 126, in import_module
    return _bootstrap._gcd_import(name[level:], package, level)
  File "<frozen importlib._bootstrap>", line 1050, in _gcd_import
  File "<frozen importlib._bootstrap>", line 1027, in _find_and_load
  File "<frozen importlib._bootstrap>", line 1006, in _find_and_load_unlocked
  File "<frozen importlib._bootstrap>", line 688, in _load_unlocked
  File "<frozen importlib._bootstrap_external>", line 883, in exec_module
  File "<frozen importlib._bootstrap>", line 241, in _call_with_frames_removed
  File "/bi/home/dkreil/.conda/envs/mercat2/lib/python3.10/site-packages/ray/dashboard/modules/dashboard_sdk.py", line 10, in <module>
    from pkg_resources import packaging
ImportError: cannot import name 'packaging' from 'pkg_resources' (/bi/home/dkreil/.conda/envs/mercat2/lib/python3.10/site-packages/pkg_resources/__init__.py)
2024-05-23 19:52:54,252 INFO worker.py:1621 -- Started a local Ray instance.
[2024-05-23 19:52:54,936 E 40696 40696] core_worker.cc:201: Failed to register worker 01000000ffffffffffffffffffffffffffffffffffffffffffffffff to Raylet. IOError: [RayletClient] Unable to register worker with raylet. No such file or directory

Error messages from mercat2 after first starting a Ray head node:

Starting MerCat2 v1.4.1 with k-mer 1 and 48 threads

2024-05-23 20:23:58,024 INFO worker.py:1431 -- Connecting to existing Ray cluster at address: 141.244.140.16:6379...
Traceback (most recent call last):
  File "/bi/home/dkreil/.conda/envs/mercat2/bin/mercat2.py", line 508, in <module>
    mercat_main()
  File "/bi/home/dkreil/.conda/envs/mercat2/bin/mercat2.py", line 217, in mercat_main
    ray.init(num_cpus=m_num_cores, log_to_driver=DEBUG)
  File "/bi/home/dkreil/.conda/envs/mercat2/lib/python3.10/site-packages/ray/_private/client_mode_hook.py", line 103, in wrapper
    return func(*args, **kwargs)
  File "/bi/home/dkreil/.conda/envs/mercat2/lib/python3.10/site-packages/ray/_private/worker.py", line 1523, in init
    raise ValueError(
ValueError: When connecting to an existing cluster, num_cpus and num_gpus must not be provided.

Can't perform the analysis

Hello, your tool seems fantastic, but I can't use it properly. 
I have some metagenomic reads samples (forward and reverse) that I have merged to have just one fasta file for each sample. 
I wanted to calculate the alfa diversity with Mercat2 to compare the samples with a kmer length of 31. 
But it seems that my workstation can't afford the computational effort. 
For some samples, it works, but when it starts with a 2GB sample, it complains that a "worker died" and the analysis fails. 
Do you have any recommendations? 

How to find out the k-mers summary of protein fasta

Hello! Mercat2 is a useful software to find out k-mers and I'm thankful for its help. I have a question that how to find out the k-mers summary of protein fasta after running Mercat2 using faa file. I just found a tsv file with a column title named 'count' and I couldn't understand the meaning of it. Hope your answer!