Repository: gbouras13/hybracter (language: python)

Issue:
[Feature request] Automatic Dnaapler fall-back or adjustable Dnaapler --evalue param
Thank you so much for this amazing tool!

**Is your feature request related to a problem? Please describe.**
Due to the nature of Snakemake, if an error occurs during any of the steps, the whole pipeline fails and the output is removed.
I find this often happens because of the Dnaapler step not being able to find an appropriate match. Unfortunately, Dnaapler's default search for DnaA is not usable for me, since (to my knowledge) [almost no archaea possess DnaA](https://www.ncbi.nlm.nih.gov/research/cog/cog/COG0593/). So I always have to use a custom Dnaapler database, where I chose RpoA, [because this is ubiquitous in all domains of life](https://www.ncbi.nlm.nih.gov/research/cog/cog/COG0202/). However, this sometimes results in not-that-specific hits, which fail the Dnaapler BLASTx cutoffs and therefore force the whole Hybracter pipeline to fail, even though the assembly is good enough for downstream usage (just not orientated).

**Describe the solution you'd like**
I could see two possible solutions:
1. make the Dnaapler --evalue param available as a Hybracter cli param: `hybracter long [...] --dnaapler-evalue 1e-9`.
2. Using the Snakemake [`try` and `attempt`](https://snakemake.readthedocs.io/en/stable/snakefiles/rules.html#dynamic-resources) options, one could implement a small helper python function that basically does "if first attempt use Dnaapler custom, ifelse second attempt use Dnaapler largest, else fail":

```python
def dnaaplerFallback(wildcards,attempt):
    if config.args.dnaapler_custom_db == "none" and attempt == 1:
        return "all"
    elif config.args.dnaapler_custom_db != "none" and attempt == 1:
        return "custom"
    elif config.args.dnaapler_custom_db != "none" and attempt == 2:
        return "largest"
    else:
        raise ValueError("Invalid attempt number or dnaapler_custom_db value")
```

Thank you for your consideration
- Richard
[Feature request] Optional automatic genome size estimation
Hi George,

I wanted to discuss the option of adding an 'auto' mode to the `-c` parameter.

In my opinion, the assembly of novel prokaryotes or at least the assembly of genomes that will be taxonomically classified AFTER they have been assembled, is a common use case of assemblers in general.
Therefore, it doesn't surprise me, that both [shovill](https://github.com/tseemann/shovill?tab=readme-ov-file#main-steps) and [dragonflye](https://github.com/rpetit3/dragonflye?tab=readme-ov-file#main-steps) contain a genome size estimation step.

Both shovill and dragonflye use [KMC](https://github.com/refresh-bio/KMC) for genome size estimation based on k-mer counting (see [Line 136 in shovill](https://github.com/tseemann/shovill/blob/f64d19e3fc3d0d29b1a5f6d20d10a08154f03305/bin/shovill#L136) and [Line 214 in dragonflye](https://github.com/rpetit3/dragonflye/blob/0e8be75b1f4a975f4fb7f8fd6f5ac07fd9dc1324/bin/dragonflye#L214)). Both use the basic but well known formula to estimate the genome size:
`For a given sequence of length L,  and a k-mer size of k, the total k-mer’s possible will be given by ( L – k ) + 1` (see [this paper](https://link.springer.com/protocol/10.1007/978-1-0716-3226-0_4#Sec2) or for example [this tutorial](https://bioinformatics.uconn.edu/genome-size-estimation-tutorial/)).

A way to implement this in snakemake would be something like this:

```bash
# Since kmc is unable to deal with tabs in the fastq headers (see https://github.com/refresh-bio/KMC/issues/42#issuecomment-340614580),
# we need to sanitize the input reads before running kmc. So here we replace tabs with underscores. Maybe hybracter does this already
rule sanitizeInputReads:
    input:
        longReads = lookup(
            query="sampleID == '{sampleID}'",
            within=sampleDF,
            cols="longReadFastq",
        ),
    output:
        longReads = temp(os.path.join(INTERIMPATH, "{sampleID}", "sanitized", "{sampleID}_sanitized.fastq.gz")),
    conda:
        os.path.join(workflow.basedir, "envs", "assembly.yaml")
    shell:
        """
        seqkit replace -p "\t" -r '_' {input.longReads} -o {output.longReads}
        """

# this is the actual rule that runs kmc. It needs to be defined as a checkpoint, so that the helper function can access the output easily
checkpoint kmc:
    input:
        longReads = rules.sanitizeInputReads.output.longReads,
    output:
        kmcLOG = os.path.join(RESULTPATH,"{sampleID}", "genomeSize", "{sampleID}_kmcLOG.txt"),
        dir = temp(directory(os.path.join(RESULTPATH,"{sampleID}", "genomeSize"))),
    shadow: "shallow"
    conda:
        os.path.join(workflow.basedir, "envs", "assembly.yaml")
    shell:
        "kmc -k21 -ci10 -v -fq {input.longReads} {wildcards.sampleID} {output.dir} > {output.kmcLOG} 2>&1"

# This is the helper function that opens the kmc log file and reads the number of unique k-mers
# For a given sequence of length L,  and a k-mer size of k, the total k-mer’s possible will be given by ( L – k ) + 1.
# This can be used to estimate the genome size by counting the number of unique k-mers in the sample.
# Since the genome length is much bigger than the k-mer size (L >> 21), the genome size can be estimated as the number of unique k-mers.
# Since the the estimated genome size is used to check if the chromosome assembly is complete or not, we use 80% of the estimated genome size to be on the safe side. This could obviously be adjusted.
def get_kmers(sampleID):
    with checkpoints.kmc.get(sampleID=sampleID).output.kmcLOG.open() as file:
        for line in file:
            if "No. of unique counted k-mers" in line:
                return int(float(re.search(r"No. of unique counted k-mers\s*:\s*([\d\.eE+-]+)", line).group(1)) * 0.8)


# Yes, I am running a snakemake/snaketool workflow within a snakemake workflow. Here I just want to use it as an example on how to use the helper function (see params: chromosomeSize)
rule hybracter_long:
    input:
        longReads=lookup(
                query="sampleID == '{sampleID}' & group == 'ONTonly'",
                within=sampleDF,
                cols="longReadFastq",
        ),  
# this is added, so that kmc is run before this rule
        kmc= lambda wildcards: checkpoints.kmc.get(sampleID=wildcards.sampleID).output.kmcLOG,
    output:
        dir=directory(os.path.join(RESULTPATH, "{sampleID}", "hybracter_long")),
        summary=os.path.join(RESULTPATH, "{sampleID}", "hybracter_long", "FINAL_OUTPUT", "hybracter_summary.tsv"),
    threads: workflow.cores * 1
    log:
        os.path.join(LOGPATH, "{sampleID}", "logs", "{sampleID}_hybracter_long.log"),
    conda:
        os.path.join(workflow.basedir, "envs", "assembly.yaml")
    params:
# this uses the helper function to get the genome size (or 80% of that number, see above)
        chromosomeSize = lambda wildcards: str(get_kmers(wildcards.sampleID)),
        plassemblerDatabase=config["tools"]["plassembler"]["dbpath"],
    shell:
        """
        hybracter long-single \
        -l {input.longReads} \
        -s {wildcards.sampleID} \
        -c {params.chromosomeSize} \
        -o {output.dir} \
        -t {threads} \
        --databases {params.plassemblerDatabase} \
        --no_medaka > {log} 2>&1
        """
```

I know and appreciate that you want to keep hybracter as simple and stable as possible, and adding an additional dependency like KMC would add complexity. However, I think that this could add a lot of value in terms of being able to automate the assembly process. Also, KMC is actively developed, so any bugs that could cause problems will also be fixed.

Best wishes,
Richard
[Feature request] Optional automatic genome size estimation
Hi George,

I wanted to discuss the option of adding an 'auto' mode to the `-c` parameter.

In my opinion, the assembly of novel prokaryotes or at least the assembly of genomes that will be taxonomically classified AFTER they have been assembled, is a common use case of assemblers in general.
Therefore, it doesn't surprise me, that both [shovill](https://github.com/tseemann/shovill?tab=readme-ov-file#main-steps) and [dragonflye](https://github.com/rpetit3/dragonflye?tab=readme-ov-file#main-steps) contain a genome size estimation step.

Both shovill and dragonflye use [KMC](https://github.com/refresh-bio/KMC) for genome size estimation based on k-mer counting (see [Line 136 in shovill](https://github.com/tseemann/shovill/blob/f64d19e3fc3d0d29b1a5f6d20d10a08154f03305/bin/shovill#L136) and [Line 214 in dragonflye](https://github.com/rpetit3/dragonflye/blob/0e8be75b1f4a975f4fb7f8fd6f5ac07fd9dc1324/bin/dragonflye#L214)). Both use the basic but well known formula to estimate the genome size:
`For a given sequence of length L,  and a k-mer size of k, the total k-mer’s possible will be given by ( L – k ) + 1` (see [this paper](https://link.springer.com/protocol/10.1007/978-1-0716-3226-0_4#Sec2) or for example [this tutorial](https://bioinformatics.uconn.edu/genome-size-estimation-tutorial/)).

A way to implement this in snakemake would be something like this:

```bash
# Since kmc is unable to deal with tabs in the fastq headers (see https://github.com/refresh-bio/KMC/issues/42#issuecomment-340614580),
# we need to sanitize the input reads before running kmc. So here we replace tabs with underscores. Maybe hybracter does this already
rule sanitizeInputReads:
    input:
        longReads = lookup(
            query="sampleID == '{sampleID}'",
            within=sampleDF,
            cols="longReadFastq",
        ),
    output:
        longReads = temp(os.path.join(INTERIMPATH, "{sampleID}", "sanitized", "{sampleID}_sanitized.fastq.gz")),
    conda:
        os.path.join(workflow.basedir, "envs", "assembly.yaml")
    shell:
        """
        seqkit replace -p "\t" -r '_' {input.longReads} -o {output.longReads}
        """

# this is the actual rule that runs kmc. It needs to be defined as a checkpoint, so that the helper function can access the output easily
checkpoint kmc:
    input:
        longReads = rules.sanitizeInputReads.output.longReads,
    output:
        kmcLOG = os.path.join(RESULTPATH,"{sampleID}", "genomeSize", "{sampleID}_kmcLOG.txt"),
        dir = temp(directory(os.path.join(RESULTPATH,"{sampleID}", "genomeSize"))),
    shadow: "shallow"
    conda:
        os.path.join(workflow.basedir, "envs", "assembly.yaml")
    shell:
        "kmc -k21 -ci10 -v -fq {input.longReads} {wildcards.sampleID} {output.dir} > {output.kmcLOG} 2>&1"

# This is the helper function that opens the kmc log file and reads the number of unique k-mers
# For a given sequence of length L,  and a k-mer size of k, the total k-mer’s possible will be given by ( L – k ) + 1.
# This can be used to estimate the genome size by counting the number of unique k-mers in the sample.
# Since the genome length is much bigger than the k-mer size (L >> 21), the genome size can be estimated as the number of unique k-mers.
# Since the the estimated genome size is used to check if the chromosome assembly is complete or not, we use 80% of the estimated genome size to be on the safe side. This could obviously be adjusted.
def get_kmers(sampleID):
    with checkpoints.kmc.get(sampleID=sampleID).output.kmcLOG.open() as file:
        for line in file:
            if "No. of unique counted k-mers" in line:
                return int(float(re.search(r"No. of unique counted k-mers\s*:\s*([\d\.eE+-]+)", line).group(1)) * 0.8)


# Yes, I am running a snakemake/snaketool workflow within a snakemake workflow. Here I just want to use it as an example on how to use the helper function (see params: chromosomeSize)
rule hybracter_long:
    input:
        longReads=lookup(
                query="sampleID == '{sampleID}' & group == 'ONTonly'",
                within=sampleDF,
                cols="longReadFastq",
        ),  
# this is added, so that kmc is run before this rule
        kmc= lambda wildcards: checkpoints.kmc.get(sampleID=wildcards.sampleID).output.kmcLOG,
    output:
        dir=directory(os.path.join(RESULTPATH, "{sampleID}", "hybracter_long")),
        summary=os.path.join(RESULTPATH, "{sampleID}", "hybracter_long", "FINAL_OUTPUT", "hybracter_summary.tsv"),
    threads: workflow.cores * 1
    log:
        os.path.join(LOGPATH, "{sampleID}", "logs", "{sampleID}_hybracter_long.log"),
    conda:
        os.path.join(workflow.basedir, "envs", "assembly.yaml")
    params:
# this uses the helper function to get the genome size (or 80% of that number, see above)
        chromosomeSize = lambda wildcards: str(get_kmers(wildcards.sampleID)),
        plassemblerDatabase=config["tools"]["plassembler"]["dbpath"],
    shell:
        """
        hybracter long-single \
        -l {input.longReads} \
        -s {wildcards.sampleID} \
        -c {params.chromosomeSize} \
        -o {output.dir} \
        -t {threads} \
        --databases {params.plassemblerDatabase} \
        --no_medaka > {log} 2>&1
        """
```

I know and appreciate that you want to keep hybracter as simple and stable as possible, and adding an additional dependency like KMC would add complexity. However, I think that this could add a lot of value in terms of being able to automate the assembly process. Also, KMC is actively developed, so any bugs that could cause problems will also be fixed.

Best wishes,
Richard


Discussion / hints:
Hi @richardstoeckl ,

Thanks for your issue. Very interesting, I haven't run Hybracter on archaeal genomes so I didn't know this!

I think I will do 2 things. 

The first is to implement your first suggestion. The reason for this is because I don't want to make hybracter more fragile and more complicated with fallbacks (I foresee the main problem/source of issues going forward is fragility!), but adding `--dnaapler_evalue` would be trivial, which a fall back mechanism may do more harm than good.

The other thing is that clearly you (and surely others) would find it useful to have a database for Archaea. So I'll add that in at some point soon and release an update :)

George

 



Hi @richardstoeckl ,

Have run hybracter since the Dnaapler upgrade to v0.8.0 on your archaea? I would be interested to see how it goes - and whether you still think it is worth that I add `--dnaapler_evalue` 

George
Hi @gbouras13 ,

sorry, I didn't catch a hybracter release with the new version of Dnaapler. Which version can I use to test? 

Thanks,
Richard
@richardstoeckl 

The existing version of hybracter v0.7.3 will install the newest dnaapler version (in the way snakemake works). So if either you reinstall hybracter v0.7.3, or reinstall ust the dnaapler env or even just install dnaapler v0.8.0 and run that manually, that should work. 

George
@gbouras13 

Everything seems fine with 90% of my samples, but two of them return the Error "return code -11":

```bash
Error in rule dnaapler_no_medaka:
    jobid: 28
    input: results/BBR080903_cB0304/hybracter_long/processing/pre_polish/BBR080903_cB0304_chromosome.fasta, results/BBR080903_cB0304/hybracter_long/processing/pre_polish/BBR080903_cB0304_ignore_list.txt
    output: results/BBR080903_cB0304/hybracter_long/processing/complete/dnaapler/BBR080903_cB0304/BBR080903_cB0304_reoriented.fasta, results/BBR080903_cB0304/hybracter_long/versions/BBR080903_cB0304/dnaapler.version
    log: results/BBR080903_cB0304/hybracter_long/stderr/dnaapler/BBR080903_cB0304.log (check log file(s) for error details)
    conda-env: /home/rks/phd/ouroboros/.snakemake/conda/37009296e852720e3812295fad89acac_/lib/python3.12/site-packages/hybracter/workflow/conda/61381f58265b1d3b2f5a23d599a4e749_
    shell:
        
        dnaapler all -i results/BBR080903_cB0304/hybracter_long/processing/pre_polish/BBR080903_cB0304_chromosome.fasta -o results/BBR080903_cB0304/hybracter_long/processing/complete/dnaapler/BBR080903_cB0304 --ignore results/BBR080903_cB0304/hybracter_long/processing/pre_polish/BBR080903_cB0304_ignore_list.txt -p BBR080903_cB0304 -t 8 -a nearest --db dnaa,repa -f 2> results/BBR080903_cB0304/hybracter_long/stderr/dnaapler/BBR080903_cB0304.log
        dnaapler --version > results/BBR080903_cB0304/hybracter_long/versions/BBR080903_cB0304/dnaapler.version
        
        (one of the commands exited with non-zero exit code; note that snakemake uses bash strict mode!)

Logfile results/BBR080903_cB0304/hybracter_long/stderr/dnaapler/BBR080903_cB0304.log:
================================================================================
2024-08-20 13:21:10.900 | INFO     | dnaapler.utils.validation:instantiate_dirs:23 - Checking the output directory results/BBR080903_cB0304/hybracter_long/processing/complete/dnaapler/BBR080903_cB0304
2024-08-20 13:21:10.908 | INFO     | dnaapler.utils.util:begin_dnaapler:71 - You are using dnaapler version 0.8.0
2024-08-20 13:21:10.908 | INFO     | dnaapler.utils.util:begin_dnaapler:72 - Repository homepage is https://github.com/gbouras13/dnaapler
2024-08-20 13:21:10.908 | INFO     | dnaapler.utils.util:begin_dnaapler:73 - Written by George Bouras: george.bouras@adelaide.edu.au
2024-08-20 13:21:10.908 | INFO     | dnaapler.utils.util:begin_dnaapler:74 - Your input FASTA is results/BBR080903_cB0304/hybracter_long/processing/pre_polish/BBR080903_cB0304_chromosome.fasta
2024-08-20 13:21:10.908 | INFO     | dnaapler.utils.util:begin_dnaapler:75 - Your output directory  is results/BBR080903_cB0304/hybracter_long/processing/complete/dnaapler/BBR080903_cB0304
2024-08-20 13:21:10.908 | INFO     | dnaapler.utils.util:begin_dnaapler:76 - You have specified 8 threads to use with blastx
2024-08-20 13:21:10.908 | INFO     | dnaapler.utils.util:begin_dnaapler:77 - You have specified dnaA,repA gene(s) to reorient your sequence
2024-08-20 13:21:10.908 | INFO     | dnaapler.utils.util:check_blast_version:117 - Checking BLAST installation.
2024-08-20 13:21:11.003 | INFO     | dnaapler.utils.util:check_blast_version:137 - BLAST version found is v2.16.0.
2024-08-20 13:21:11.004 | INFO     | dnaapler.utils.util:check_blast_version:147 - BLAST version is ok.
2024-08-20 13:21:11.004 | INFO     | dnaapler.utils.util:check_pyrodigal_version:92 - Checking pyrodigal installation.
2024-08-20 13:21:11.004 | INFO     | dnaapler.utils.util:check_pyrodigal_version:103 - Pyrodigal version is v3.5.1
2024-08-20 13:21:11.005 | INFO     | dnaapler.utils.util:check_pyrodigal_version:104 - Pyrodigal version is ok.
2024-08-20 13:21:11.005 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --input results/BBR080903_cB0304/hybracter_long/processing/pre_polish/BBR080903_cB0304_chromosome.fasta
2024-08-20 13:21:11.005 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --output results/BBR080903_cB0304/hybracter_long/processing/complete/dnaapler/BBR080903_cB0304
2024-08-20 13:21:11.006 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --threads 8
2024-08-20 13:21:11.006 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --force True
2024-08-20 13:21:11.006 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --prefix BBR080903_cB0304
2024-08-20 13:21:11.006 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --evalue 1e-10
2024-08-20 13:21:11.007 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --autocomplete nearest
2024-08-20 13:21:11.007 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --seed_value 13
2024-08-20 13:21:11.007 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --ignore results/BBR080903_cB0304/hybracter_long/processing/pre_polish/BBR080903_cB0304_ignore_list.txt
2024-08-20 13:21:11.007 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --custom_db 
2024-08-20 13:21:11.008 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --db dnaa,repa
2024-08-20 13:21:11.008 | INFO     | dnaapler.utils.validation:validate_fasta_all:100 - Checking that the input file results/BBR080903_cB0304/hybracter_long/processing/pre_polish/BBR080903_cB0304_chromosome.fasta is in FASTA format and has at least 1 entry.
2024-08-20 13:21:11.033 | INFO     | dnaapler.utils.validation:validate_fasta_all:107 - results/BBR080903_cB0304/hybracter_long/processing/pre_polish/BBR080903_cB0304_chromosome.fasta file checked.
2024-08-20 13:21:11.041 | INFO     | dnaapler.utils.validation:validate_fasta_all:116 - results/BBR080903_cB0304/hybracter_long/processing/pre_polish/BBR080903_cB0304_chromosome.fasta has only one entry.
2024-08-20 13:21:11.046 | INFO     | dnaapler.utils.validation:check_evalue:187 - You have specified an evalue of 1e-10.
2024-08-20 13:21:11.047 | INFO     | dnaapler:all:930 - You have specified contigs to ignore in results/BBR080903_cB0304/hybracter_long/processing/pre_polish/BBR080903_cB0304_ignore_list.txt.
2024-08-20 13:21:11.047 | INFO     | dnaapler.utils.external_tools:run:53 - Started running blastx -db /home/rks/phd/ouroboros/.snakemake/conda/37009296e852720e3812295fad89acac_/lib/python3.12/site-packages/hybracter/workflow/conda/61381f58265b1d3b2f5a23d599a4e749_/lib/python3.12/site-packages/dnaapler/db/dnaA_repA_db -evalue 1e-10 -num_threads 8 -outfmt ' 6 qseqid qlen sseqid slen length qstart qend sstart send pident nident gaps mismatch evalue bitscore qseq sseq ' -out results/BBR080903_cB0304/hybracter_long/processing/complete/dnaapler/BBR080903_cB0304/BBR080903_cB0304_blast_output.txt -query results/BBR080903_cB0304/hybracter_long/processing/pre_polish/BBR080903_cB0304_chromosome.fasta ...
2024-08-20 13:21:27.365 | ERROR    | dnaapler.utils.external_tools:run_tool:91 - Error calling blastx -db /home/rks/phd/ouroboros/.snakemake/conda/37009296e852720e3812295fad89acac_/lib/python3.12/site-packages/hybracter/workflow/conda/61381f58265b1d3b2f5a23d599a4e749_/lib/python3.12/site-packages/dnaapler/db/dnaA_repA_db -evalue 1e-10 -num_threads 8 -outfmt ' 6 qseqid qlen sseqid slen length qstart qend sstart send pident nident gaps mismatch evalue bitscore qseq sseq ' -out results/BBR080903_cB0304/hybracter_long/processing/complete/dnaapler/BBR080903_cB0304/BBR080903_cB0304_blast_output.txt -query results/BBR080903_cB0304/hybracter_long/processing/pre_polish/BBR080903_cB0304_chromosome.fasta (return code -11)
================================================================================
```
and

```bash
[Mon Aug 19 19:26:17 2024]
Error in rule dnaapler_no_medaka:
    jobid: 28
    input: results/Shy_cB0304/hybracter_long/processing/pre_polish/Shy_cB0304_chromosome.fasta, results/Shy_cB0304/hybracter_long/processing/pre_polish/Shy_cB0304_ignore_list.txt
    output: results/Shy_cB0304/hybracter_long/processing/complete/dnaapler/Shy_cB0304/Shy_cB0304_reoriented.fasta, results/Shy_cB0304/hybracter_long/versions/Shy_cB0304/dnaapler.version
    log: results/Shy_cB0304/hybracter_long/stderr/dnaapler/Shy_cB0304.log (check log file(s) for error details)
    conda-env: /home/rks/phd/ouroboros/.snakemake/conda/37009296e852720e3812295fad89acac_/lib/python3.12/site-packages/hybracter/workflow/conda/61381f58265b1d3b2f5a23d599a4e749_
    shell:
        
        dnaapler all -i results/Shy_cB0304/hybracter_long/processing/pre_polish/Shy_cB0304_chromosome.fasta -o results/Shy_cB0304/hybracter_long/processing/complete/dnaapler/Shy_cB0304 --ignore results/Shy_cB0304/hybracter_long/processing/pre_polish/Shy_cB0304_ignore_list.txt -p Shy_cB0304 -t 8 -a nearest --db dnaa,repa -f 2> results/Shy_cB0304/hybracter_long/stderr/dnaapler/Shy_cB0304.log
        dnaapler --version > results/Shy_cB0304/hybracter_long/versions/Shy_cB0304/dnaapler.version
        
        (one of the commands exited with non-zero exit code; note that snakemake uses bash strict mode!)
Logfile results/Shy_cB0304/hybracter_long/stderr/dnaapler/Shy_cB0304.log:
================================================================================
2024-08-19 19:26:04.054 | INFO     | dnaapler.utils.validation:instantiate_dirs:23 - Checking the output directory results/Shy_cB0304/hybracter_long/processing/complete/dnaapler/Shy_cB0304
2024-08-19 19:26:04.061 | INFO     | dnaapler.utils.util:begin_dnaapler:71 - You are using dnaapler version 0.8.0
2024-08-19 19:26:04.062 | INFO     | dnaapler.utils.util:begin_dnaapler:72 - Repository homepage is https://github.com/gbouras13/dnaapler
2024-08-19 19:26:04.062 | INFO     | dnaapler.utils.util:begin_dnaapler:73 - Written by George Bouras: george.bouras@adelaide.edu.au
2024-08-19 19:26:04.062 | INFO     | dnaapler.utils.util:begin_dnaapler:74 - Your input FASTA is results/Shy_cB0304/hybracter_long/processing/pre_polish/Shy_cB0304_chromosome.fasta
2024-08-19 19:26:04.062 | INFO     | dnaapler.utils.util:begin_dnaapler:75 - Your output directory  is results/Shy_cB0304/hybracter_long/processing/complete/dnaapler/Shy_cB0304
2024-08-19 19:26:04.062 | INFO     | dnaapler.utils.util:begin_dnaapler:76 - You have specified 8 threads to use with blastx
2024-08-19 19:26:04.062 | INFO     | dnaapler.utils.util:begin_dnaapler:77 - You have specified dnaA,repA gene(s) to reorient your sequence
2024-08-19 19:26:04.062 | INFO     | dnaapler.utils.util:check_blast_version:117 - Checking BLAST installation.
2024-08-19 19:26:04.147 | INFO     | dnaapler.utils.util:check_blast_version:137 - BLAST version found is v2.16.0.
2024-08-19 19:26:04.148 | INFO     | dnaapler.utils.util:check_blast_version:147 - BLAST version is ok.
2024-08-19 19:26:04.148 | INFO     | dnaapler.utils.util:check_pyrodigal_version:92 - Checking pyrodigal installation.
2024-08-19 19:26:04.148 | INFO     | dnaapler.utils.util:check_pyrodigal_version:103 - Pyrodigal version is v3.5.1
2024-08-19 19:26:04.148 | INFO     | dnaapler.utils.util:check_pyrodigal_version:104 - Pyrodigal version is ok.
2024-08-19 19:26:04.148 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --input results/Shy_cB0304/hybracter_long/processing/pre_polish/Shy_cB0304_chromosome.fasta
2024-08-19 19:26:04.149 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --output results/Shy_cB0304/hybracter_long/processing/complete/dnaapler/Shy_cB0304
2024-08-19 19:26:04.149 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --threads 8
2024-08-19 19:26:04.149 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --force True
2024-08-19 19:26:04.149 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --prefix Shy_cB0304
2024-08-19 19:26:04.149 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --evalue 1e-10
2024-08-19 19:26:04.149 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --autocomplete nearest
2024-08-19 19:26:04.149 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --seed_value 13
2024-08-19 19:26:04.149 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --ignore results/Shy_cB0304/hybracter_long/processing/pre_polish/Shy_cB0304_ignore_list.txt
2024-08-19 19:26:04.149 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --custom_db 
2024-08-19 19:26:04.150 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --db dnaa,repa
2024-08-19 19:26:04.150 | INFO     | dnaapler.utils.validation:validate_fasta_all:100 - Checking that the input file results/Shy_cB0304/hybracter_long/processing/pre_polish/Shy_cB0304_chromosome.fasta is in FASTA format and has at least 1 entry.
2024-08-19 19:26:04.162 | INFO     | dnaapler.utils.validation:validate_fasta_all:107 - results/Shy_cB0304/hybracter_long/processing/pre_polish/Shy_cB0304_chromosome.fasta file checked.
2024-08-19 19:26:04.168 | INFO     | dnaapler.utils.validation:validate_fasta_all:116 - results/Shy_cB0304/hybracter_long/processing/pre_polish/Shy_cB0304_chromosome.fasta has only one entry.
2024-08-19 19:26:04.173 | INFO     | dnaapler.utils.validation:check_evalue:187 - You have specified an evalue of 1e-10.
2024-08-19 19:26:04.173 | INFO     | dnaapler:all:930 - You have specified contigs to ignore in results/Shy_cB0304/hybracter_long/processing/pre_polish/Shy_cB0304_ignore_list.txt.
2024-08-19 19:26:04.174 | INFO     | dnaapler.utils.external_tools:run:53 - Started running blastx -db /home/rks/phd/ouroboros/.snakemake/conda/37009296e852720e3812295fad89acac_/lib/python3.12/site-packages/hybracter/workflow/conda/61381f58265b1d3b2f5a23d599a4e749_/lib/python3.12/site-packages/dnaapler/db/dnaA_repA_db -evalue 1e-10 -num_threads 8 -outfmt ' 6 qseqid qlen sseqid slen length qstart qend sstart send pident nident gaps mismatch evalue bitscore qseq sseq ' -out results/Shy_cB0304/hybracter_long/processing/complete/dnaapler/Shy_cB0304/Shy_cB0304_blast_output.txt -query results/Shy_cB0304/hybracter_long/processing/pre_polish/Shy_cB0304_chromosome.fasta ...
2024-08-19 19:26:16.998 | ERROR    | dnaapler.utils.external_tools:run_tool:91 - Error calling blastx -db /home/rks/phd/ouroboros/.snakemake/conda/37009296e852720e3812295fad89acac_/lib/python3.12/site-packages/hybracter/workflow/conda/61381f58265b1d3b2f5a23d599a4e749_/lib/python3.12/site-packages/dnaapler/db/dnaA_repA_db -evalue 1e-10 -num_threads 8 -outfmt ' 6 qseqid qlen sseqid slen length qstart qend sstart send pident nident gaps mismatch evalue bitscore qseq sseq ' -out results/Shy_cB0304/hybracter_long/processing/complete/dnaapler/Shy_cB0304/Shy_cB0304_blast_output.txt -query results/Shy_cB0304/hybracter_long/processing/pre_polish/Shy_cB0304_chromosome.fasta (return code -11)
================================================================================

[Mon Aug 19 19:26:17 2024]
Error in rule dnaapler_pre_chrom:
    jobid: 33
    input: results/Shy_cB0304/hybracter_long/processing/pre_polish/Shy_cB0304_chromosome.fasta, results/Shy_cB0304/hybracter_long/processing/pre_polish/Shy_cB0304_ignore_list.txt
    output: results/Shy_cB0304/hybracter_long/processing/complete/dnaapler/Shy_cB0304_pre_chrom/Shy_cB0304_reoriented.fasta
    log: results/Shy_cB0304/hybracter_long/stderr/dnaapler/Shy_cB0304_pre_chrom.log (check log file(s) for error details)
    conda-env: /home/rks/phd/ouroboros/.snakemake/conda/37009296e852720e3812295fad89acac_/lib/python3.12/site-packages/hybracter/workflow/conda/61381f58265b1d3b2f5a23d599a4e749_
    shell:
        
        dnaapler all -i results/Shy_cB0304/hybracter_long/processing/pre_polish/Shy_cB0304_chromosome.fasta -o results/Shy_cB0304/hybracter_long/processing/complete/dnaapler/Shy_cB0304_pre_chrom --ignore results/Shy_cB0304/hybracter_long/processing/pre_polish/Shy_cB0304_ignore_list.txt -p Shy_cB0304 -t 8 -a nearest --db dnaa,repa -f 2> results/Shy_cB0304/hybracter_long/stderr/dnaapler/Shy_cB0304_pre_chrom.log
        
        (one of the commands exited with non-zero exit code; note that snakemake uses bash strict mode!)
Logfile results/Shy_cB0304/hybracter_long/stderr/dnaapler/Shy_cB0304_pre_chrom.log:
================================================================================
2024-08-19 19:26:04.063 | INFO     | dnaapler.utils.validation:instantiate_dirs:23 - Checking the output directory results/Shy_cB0304/hybracter_long/processing/complete/dnaapler/Shy_cB0304_pre_chrom
2024-08-19 19:26:04.071 | INFO     | dnaapler.utils.util:begin_dnaapler:71 - You are using dnaapler version 0.8.0
2024-08-19 19:26:04.071 | INFO     | dnaapler.utils.util:begin_dnaapler:72 - Repository homepage is https://github.com/gbouras13/dnaapler
2024-08-19 19:26:04.071 | INFO     | dnaapler.utils.util:begin_dnaapler:73 - Written by George Bouras: george.bouras@adelaide.edu.au
2024-08-19 19:26:04.071 | INFO     | dnaapler.utils.util:begin_dnaapler:74 - Your input FASTA is results/Shy_cB0304/hybracter_long/processing/pre_polish/Shy_cB0304_chromosome.fasta
2024-08-19 19:26:04.071 | INFO     | dnaapler.utils.util:begin_dnaapler:75 - Your output directory  is results/Shy_cB0304/hybracter_long/processing/complete/dnaapler/Shy_cB0304_pre_chrom
2024-08-19 19:26:04.071 | INFO     | dnaapler.utils.util:begin_dnaapler:76 - You have specified 8 threads to use with blastx
2024-08-19 19:26:04.071 | INFO     | dnaapler.utils.util:begin_dnaapler:77 - You have specified dnaA,repA gene(s) to reorient your sequence
2024-08-19 19:26:04.071 | INFO     | dnaapler.utils.util:check_blast_version:117 - Checking BLAST installation.
2024-08-19 19:26:04.148 | INFO     | dnaapler.utils.util:check_blast_version:137 - BLAST version found is v2.16.0.
2024-08-19 19:26:04.148 | INFO     | dnaapler.utils.util:check_blast_version:147 - BLAST version is ok.
2024-08-19 19:26:04.149 | INFO     | dnaapler.utils.util:check_pyrodigal_version:92 - Checking pyrodigal installation.
2024-08-19 19:26:04.149 | INFO     | dnaapler.utils.util:check_pyrodigal_version:103 - Pyrodigal version is v3.5.1
2024-08-19 19:26:04.149 | INFO     | dnaapler.utils.util:check_pyrodigal_version:104 - Pyrodigal version is ok.
2024-08-19 19:26:04.149 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --input results/Shy_cB0304/hybracter_long/processing/pre_polish/Shy_cB0304_chromosome.fasta
2024-08-19 19:26:04.150 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --output results/Shy_cB0304/hybracter_long/processing/complete/dnaapler/Shy_cB0304_pre_chrom
2024-08-19 19:26:04.150 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --threads 8
2024-08-19 19:26:04.150 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --force True
2024-08-19 19:26:04.151 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --prefix Shy_cB0304
2024-08-19 19:26:04.151 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --evalue 1e-10
2024-08-19 19:26:04.151 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --autocomplete nearest
2024-08-19 19:26:04.151 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --seed_value 13
2024-08-19 19:26:04.151 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --ignore results/Shy_cB0304/hybracter_long/processing/pre_polish/Shy_cB0304_ignore_list.txt
2024-08-19 19:26:04.151 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --custom_db 
2024-08-19 19:26:04.151 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --db dnaa,repa
2024-08-19 19:26:04.151 | INFO     | dnaapler.utils.validation:validate_fasta_all:100 - Checking that the input file results/Shy_cB0304/hybracter_long/processing/pre_polish/Shy_cB0304_chromosome.fasta is in FASTA format and has at least 1 entry.
2024-08-19 19:26:04.163 | INFO     | dnaapler.utils.validation:validate_fasta_all:107 - results/Shy_cB0304/hybracter_long/processing/pre_polish/Shy_cB0304_chromosome.fasta file checked.
2024-08-19 19:26:04.169 | INFO     | dnaapler.utils.validation:validate_fasta_all:116 - results/Shy_cB0304/hybracter_long/processing/pre_polish/Shy_cB0304_chromosome.fasta has only one entry.
2024-08-19 19:26:04.174 | INFO     | dnaapler.utils.validation:check_evalue:187 - You have specified an evalue of 1e-10.
2024-08-19 19:26:04.174 | INFO     | dnaapler:all:930 - You have specified contigs to ignore in results/Shy_cB0304/hybracter_long/processing/pre_polish/Shy_cB0304_ignore_list.txt.
2024-08-19 19:26:04.175 | INFO     | dnaapler.utils.external_tools:run:53 - Started running blastx -db /home/rks/phd/ouroboros/.snakemake/conda/37009296e852720e3812295fad89acac_/lib/python3.12/site-packages/hybracter/workflow/conda/61381f58265b1d3b2f5a23d599a4e749_/lib/python3.12/site-packages/dnaapler/db/dnaA_repA_db -evalue 1e-10 -num_threads 8 -outfmt ' 6 qseqid qlen sseqid slen length qstart qend sstart send pident nident gaps mismatch evalue bitscore qseq sseq ' -out results/Shy_cB0304/hybracter_long/processing/complete/dnaapler/Shy_cB0304_pre_chrom/Shy_cB0304_blast_output.txt -query results/Shy_cB0304/hybracter_long/processing/pre_polish/Shy_cB0304_chromosome.fasta ...
2024-08-19 19:26:17.117 | ERROR    | dnaapler.utils.external_tools:run_tool:91 - Error calling blastx -db /home/rks/phd/ouroboros/.snakemake/conda/37009296e852720e3812295fad89acac_/lib/python3.12/site-packages/hybracter/workflow/conda/61381f58265b1d3b2f5a23d599a4e749_/lib/python3.12/site-packages/dnaapler/db/dnaA_repA_db -evalue 1e-10 -num_threads 8 -outfmt ' 6 qseqid qlen sseqid slen length qstart qend sstart send pident nident gaps mismatch evalue bitscore qseq sseq ' -out results/Shy_cB0304/hybracter_long/processing/complete/dnaapler/Shy_cB0304_pre_chrom/Shy_cB0304_blast_output.txt -query results/Shy_cB0304/hybracter_long/processing/pre_polish/Shy_cB0304_chromosome.fasta (return code -11)
================================================================================
```

I will look into these samples in more detail and report back to you. Could you tell me if `dnaapler.utils.util:begin_dnaapler:82 - Parameter: --db dnaa,repa` is the intended option (so no mention of the archaeal DB), and what the blastx `(return code -11)` could refer to?
Hi @richardstoeckl ,

You are right - the dnaa,repa needs to be changed to include the COG1474 gene. My bad.

The return code is very strange, not sure about that (my guess it is related to finding nothing). 

Do you mind trying `dnaapler all` (by itself) on these contigs? Specifically the contig `results/Shy_cB0304/hybracter_long/processing/pre_polish/Shy_cB0304_chromosome.fasta`. It now includes archaea support.

I might even introduce an `--archaea` flag into hybracter that only tries COG1474. 

The reason is is that if all your samples bar 2 worked, that means a dnaA homolog was found in those. Which is quite high.

Therefore if we want the standard going forward to be the COG1474 for archaea, then we should probably make that the only option (i.e. to remove dnaa).

George






Hi @gbouras13 ,

please excuse the long time it took me to get back to this issue.

I think we now figured out and fixed the first problem:

Using the **archae mode** on the "Shy" sample worked:
```bash
dnaapler archaea -i Shy_cB0304_manual/processing/pre_polish/Shy_cB0304_manual_chromosome.fasta -o Shy_cB0304_manual/manualDnaapler/ -p Shy_cB0304_manualDnaapler -t 8 -a nearest -f 2> Shy_cB0304_manual/manualDnaapler.log
```
With the following log:

```bash
2024-09-10 15:13:50.628 | INFO     | dnaapler.utils.util:begin_dnaapler:71 - You are using dnaapler version 0.8.0
2024-09-10 15:13:50.628 | INFO     | dnaapler.utils.util:begin_dnaapler:72 - Repository homepage is https://github.com/gbouras13/dnaapler
2024-09-10 15:13:50.628 | INFO     | dnaapler.utils.util:begin_dnaapler:73 - Written by George Bouras: george.bouras@adelaide.edu.au
2024-09-10 15:13:50.628 | INFO     | dnaapler.utils.util:begin_dnaapler:74 - Your input FASTA is Shy_cB0304_manual/processing/pre_polish/Shy_cB0304_manual_chromosome.fasta
2024-09-10 15:13:50.628 | INFO     | dnaapler.utils.util:begin_dnaapler:75 - Your output directory  is Shy_cB0304_manual/manualDnaapler/
2024-09-10 15:13:50.628 | INFO     | dnaapler.utils.util:begin_dnaapler:76 - You have specified 8 threads to use with blastx
2024-09-10 15:13:50.628 | INFO     | dnaapler.utils.util:begin_dnaapler:77 - You have specified cog1474 gene(s) to reorient your sequence
2024-09-10 15:13:50.628 | INFO     | dnaapler.utils.util:check_blast_version:117 - Checking BLAST installation.
2024-09-10 15:13:50.707 | INFO     | dnaapler.utils.util:check_blast_version:137 - BLAST version found is v2.16.0.
2024-09-10 15:13:50.707 | INFO     | dnaapler.utils.util:check_blast_version:147 - BLAST version is ok.
2024-09-10 15:13:50.707 | INFO     | dnaapler.utils.util:check_pyrodigal_version:92 - Checking pyrodigal installation.
2024-09-10 15:13:50.707 | INFO     | dnaapler.utils.util:check_pyrodigal_version:103 - Pyrodigal version is v3.5.1
2024-09-10 15:13:50.707 | INFO     | dnaapler.utils.util:check_pyrodigal_version:104 - Pyrodigal version is ok.
2024-09-10 15:13:50.707 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --input Shy_cB0304_manual/processing/pre_polish/Shy_cB0304_manual_chromosome.fasta
2024-09-10 15:13:50.707 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --output Shy_cB0304_manual/manualDnaapler/
2024-09-10 15:13:50.708 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --threads 8
2024-09-10 15:13:50.708 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --force True
2024-09-10 15:13:50.708 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --prefix Shy_cB0304_manualDnaapler
2024-09-10 15:13:50.708 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --evalue 1e-10
2024-09-10 15:13:50.708 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --autocomplete nearest
2024-09-10 15:13:50.708 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --seed_value 13
2024-09-10 15:13:50.708 | INFO     | dnaapler:archaea:259 - You have chosen nearest method to reorient your sequence if the BLAST based method fails.
2024-09-10 15:13:50.708 | INFO     | dnaapler.utils.validation:validate_fasta:46 - Checking that the input file Shy_cB0304_manual/processing/pre_polish/Shy_cB0304_manual_chromosome.fasta is in FASTA format and has only 1 entry.
2024-09-10 15:13:50.718 | INFO     | dnaapler.utils.validation:validate_fasta:53 - Shy_cB0304_manual/processing/pre_polish/Shy_cB0304_manual_chromosome.fasta file checked.
2024-09-10 15:13:50.723 | INFO     | dnaapler.utils.validation:validate_fasta:62 - Shy_cB0304_manual/processing/pre_polish/Shy_cB0304_manual_chromosome.fasta has only one entry.
2024-09-10 15:13:50.724 | INFO     | dnaapler.utils.validation:check_evalue:187 - You have specified an evalue of 1e-10.
2024-09-10 15:13:50.724 | INFO     | dnaapler.utils.external_tools:run:53 - Started running blastx -db /home/rks/phd/ouroboros/.snakemake/conda/37009296e852720e3812295fad89acac_/lib/python3.12/site-packages/hybracter/workflow/conda/61381f58265b1d3b2f5a23d599a4e749_/lib/python3.12/site-packages/dnaapler/db/cog1474_db -evalue 1e-10 -num_threads 8 -outfmt ' 6 qseqid qlen sseqid slen length qstart qend sstart send pident nident gaps mismatch evalue bitscore qseq sseq ' -out Shy_cB0304_manual/manualDnaapler/Shy_cB0304_manualDnaapler_blast_output.txt -query Shy_cB0304_manual/processing/pre_polish/Shy_cB0304_manual_chromosome.fasta ...
2024-09-10 15:13:57.956 | INFO     | dnaapler.utils.external_tools:run:55 - Done running blastx -db /home/rks/phd/ouroboros/.snakemake/conda/37009296e852720e3812295fad89acac_/lib/python3.12/site-packages/hybracter/workflow/conda/61381f58265b1d3b2f5a23d599a4e749_/lib/python3.12/site-packages/dnaapler/db/cog1474_db -evalue 1e-10 -num_threads 8 -outfmt ' 6 qseqid qlen sseqid slen length qstart qend sstart send pident nident gaps mismatch evalue bitscore qseq sseq ' -out Shy_cB0304_manual/manualDnaapler/Shy_cB0304_manualDnaapler_blast_output.txt -query Shy_cB0304_manual/processing/pre_polish/Shy_cB0304_manual_chromosome.fasta
2024-09-10 15:13:57.963 | INFO     | dnaapler.utils.processing:reorient_sequence:145 - cog1474 gene identified. It starts at coordinate 574403 on the forward strand in your input file.
2024-09-10 15:13:57.964 | INFO     | dnaapler.utils.processing:reorient_sequence:148 - The best hit with a valid start codon in the database is cog1474_gi|347522555|ref|YP_004780125.1|, which has length of 403 AAs.
2024-09-10 15:13:57.964 | INFO     | dnaapler.utils.processing:reorient_sequence:151 - 396 AAs were covered by the best hit, with an overall coverage of 98.26%.
2024-09-10 15:13:57.964 | INFO     | dnaapler.utils.processing:reorient_sequence:154 - 261 AAs were identical, with an overall identity of 65.91%.
2024-09-10 15:13:57.964 | INFO     | dnaapler.utils.processing:reorient_sequence:157 - Re-orienting.
2024-09-10 15:13:57.974 | INFO     | dnaapler.utils.util:end_dnaapler:160 - dnaapler has finished
2024-09-10 15:13:57.974 | INFO     | dnaapler.utils.util:end_dnaapler:161 - Elapsed time: 7.35 seconds
```

Using the **all** mode on this sample also correctly identified the correct gene and everything worked.
```bash
dnaapler all -i Shy_cB0304_manual/processing/pre_polish/Shy_cB0304_manual_chromosome.fasta -o Shy_cB0304_manual/manualDnaapler/ --ignore Shy_cB0304_manual/processing/pre_polish/Shy_cB0304_manual_ignore_list.txt -p Shy_cB0304_manualDnaapler -t 8 -a nearest -f 2> Shy_cB0304_manual/manualDnaapler.log
```
Which resulted in:
```bash
Contig	Gene_Reoriented	Start	Strand	Top_Hit	Top_Hit_Length	Covered_Length	Coverage	Identical_AAs	Identity_Percentage
contig_1	cog1474	574403	forward	cog1474_gi|347522555|ref|YP_004780125.1|	403	396	98.26	261	65.91
```

With this logfile:
```bash
2024-09-10 15:05:38.421 | INFO     | dnaapler.utils.validation:instantiate_dirs:23 - Checking the output directory Shy_cB0304_manual/manualDnaapler/
2024-09-10 15:05:38.429 | INFO     | dnaapler.utils.util:begin_dnaapler:71 - You are using dnaapler version 0.8.0
2024-09-10 15:05:38.429 | INFO     | dnaapler.utils.util:begin_dnaapler:72 - Repository homepage is https://github.com/gbouras13/dnaapler
2024-09-10 15:05:38.429 | INFO     | dnaapler.utils.util:begin_dnaapler:73 - Written by George Bouras: george.bouras@adelaide.edu.au
2024-09-10 15:05:38.429 | INFO     | dnaapler.utils.util:begin_dnaapler:74 - Your input FASTA is Shy_cB0304_manual/processing/pre_polish/Shy_cB0304_manual_chromosome.fasta
2024-09-10 15:05:38.429 | INFO     | dnaapler.utils.util:begin_dnaapler:75 - Your output directory  is Shy_cB0304_manual/manualDnaapler/
2024-09-10 15:05:38.429 | INFO     | dnaapler.utils.util:begin_dnaapler:76 - You have specified 8 threads to use with blastx
2024-09-10 15:05:38.430 | INFO     | dnaapler.utils.util:begin_dnaapler:77 - You have specified all gene(s) to reorient your sequence
2024-09-10 15:05:38.430 | INFO     | dnaapler.utils.util:check_blast_version:117 - Checking BLAST installation.
2024-09-10 15:05:38.521 | INFO     | dnaapler.utils.util:check_blast_version:137 - BLAST version found is v2.16.0.
2024-09-10 15:05:38.521 | INFO     | dnaapler.utils.util:check_blast_version:147 - BLAST version is ok.
2024-09-10 15:05:38.522 | INFO     | dnaapler.utils.util:check_pyrodigal_version:92 - Checking pyrodigal installation.
2024-09-10 15:05:38.522 | INFO     | dnaapler.utils.util:check_pyrodigal_version:103 - Pyrodigal version is v3.5.1
2024-09-10 15:05:38.522 | INFO     | dnaapler.utils.util:check_pyrodigal_version:104 - Pyrodigal version is ok.
2024-09-10 15:05:38.522 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --input Shy_cB0304_manual/processing/pre_polish/Shy_cB0304_manual_chromosome.fasta
2024-09-10 15:05:38.522 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --output Shy_cB0304_manual/manualDnaapler/
2024-09-10 15:05:38.522 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --threads 8
2024-09-10 15:05:38.523 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --force True
2024-09-10 15:05:38.523 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --prefix Shy_cB0304_manualDnaapler
2024-09-10 15:05:38.523 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --evalue 1e-10
2024-09-10 15:05:38.523 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --autocomplete nearest
2024-09-10 15:05:38.523 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --seed_value 13
2024-09-10 15:05:38.523 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --ignore Shy_cB0304_manual/processing/pre_polish/Shy_cB0304_manual_ignore_list.txt
2024-09-10 15:05:38.523 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --custom_db 
2024-09-10 15:05:38.524 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --db all
2024-09-10 15:05:38.524 | INFO     | dnaapler.utils.validation:validate_fasta_all:100 - Checking that the input file Shy_cB0304_manual/processing/pre_polish/Shy_cB0304_manual_chromosome.fasta is in FASTA format and has at least 1 entry.
2024-09-10 15:05:38.542 | INFO     | dnaapler.utils.validation:validate_fasta_all:107 - Shy_cB0304_manual/processing/pre_polish/Shy_cB0304_manual_chromosome.fasta file checked.
2024-09-10 15:05:38.549 | INFO     | dnaapler.utils.validation:validate_fasta_all:116 - Shy_cB0304_manual/processing/pre_polish/Shy_cB0304_manual_chromosome.fasta has only one entry.
2024-09-10 15:05:38.553 | INFO     | dnaapler.utils.validation:check_evalue:187 - You have specified an evalue of 1e-10.
2024-09-10 15:05:38.554 | INFO     | dnaapler:all:930 - You have specified contigs to ignore in Shy_cB0304_manual/processing/pre_polish/Shy_cB0304_manual_ignore_list.txt.
2024-09-10 15:05:38.554 | INFO     | dnaapler.utils.external_tools:run:53 - Started running blastx -db /home/rks/phd/ouroboros/.snakemake/conda/37009296e852720e3812295fad89acac_/lib/python3.12/site-packages/hybracter/workflow/conda/61381f58265b1d3b2f5a23d599a4e749_/lib/python3.12/site-packages/dnaapler/db/all_db -evalue 1e-10 -num_threads 8 -outfmt ' 6 qseqid qlen sseqid slen length qstart qend sstart send pident nident gaps mismatch evalue bitscore qseq sseq ' -out Shy_cB0304_manual/manualDnaapler/Shy_cB0304_manualDnaapler_blast_output.txt -query Shy_cB0304_manual/processing/pre_polish/Shy_cB0304_manual_chromosome.fasta ...
2024-09-10 15:06:47.829 | INFO     | dnaapler.utils.external_tools:run:55 - Done running blastx -db /home/rks/phd/ouroboros/.snakemake/conda/37009296e852720e3812295fad89acac_/lib/python3.12/site-packages/hybracter/workflow/conda/61381f58265b1d3b2f5a23d599a4e749_/lib/python3.12/site-packages/dnaapler/db/all_db -evalue 1e-10 -num_threads 8 -outfmt ' 6 qseqid qlen sseqid slen length qstart qend sstart send pident nident gaps mismatch evalue bitscore qseq sseq ' -out Shy_cB0304_manual/manualDnaapler/Shy_cB0304_manualDnaapler_blast_output.txt -query Shy_cB0304_manual/processing/pre_polish/Shy_cB0304_manual_chromosome.fasta
2024-09-10 15:06:47.829 | WARNING  | dnaapler:all:973 - Shy_cB0304_manual/processing/pre_polish/Shy_cB0304_manual_ignore_list.txt contains no text. No contigs will be ignored
2024-09-10 15:06:47.851 | INFO     | dnaapler.utils.util:end_dnaapler:160 - dnaapler has finished
2024-09-10 15:06:47.851 | INFO     | dnaapler.utils.util:end_dnaapler:161 - Elapsed time: 69.43 seconds
```


However, I think you still need to update the dnaapler command used within hybracter to remove or adjust the `--db dnaa,repa` parameter.

I will look into the other samples that DIDN'T fail more closely and get back to you!


So just looking into the very first example, I think adjusting the `--db` parameter would fix the problem.

This rule was run by using hybracter v0.7.3 (from the hybracter.log):

```bash
[2024:08:19 16:03:12] Copying system default config to results/BBR060204_cB0304/hybracter_long/config.yaml
[2024:08:19 16:03:12] Updating config file with new values
[2024:08:19 16:03:12] Writing config file to results/BBR060204_cB0304/hybracter_long/config.yaml
[2024:08:19 16:03:12] ------------------
[2024:08:19 16:03:12] | Runtime config |
[2024:08:19 16:03:12] ------------------

args:
  chromosome: 1449461
  configfile: results/BBR060204_cB0304/hybracter_long/config.yaml
  contaminants: none
  databases: /mnt/DATA/common/databases/plassembler/
  depth_filter: 0.25
  dnaapler_custom_db: none
  flyeModel: --nano-hq
  log: results/BBR060204_cB0304/hybracter_long/hybracter.log
  logic: best
  longreads: /mnt/DATA2/rks/BBRseq03AndBBRseq04combined/BBR060204_combinedBBRseq03And04.fastq.gz
  medakaModel: r1041_e82_400bps_sup_v4.2.0
  min_length: 1000
  min_quality: 9
  no_medaka: true
  output: results/BBR060204_cB0304/hybracter_long
  sample: BBR060204_cB0304
  single: true
  skip_qc: false
  subsample_depth: 100
qc:
  compression: 5
  hostRemoveFlagstat: -f 4 -F 3584
  minimapIndex: -I 8G
  minimapModel: map-ont
resources:
  big:
    cpu: 16
    mem: 32000
    time: '23:59:00'
  med:
    cpu: 8
    mem: 16000
    time: 08:00:00
  sml:
    cpu: 1
    mem: 4000
    time: 00:00:05



localrule dnaapler_no_medaka:
    input: results/BBR060204_cB0304/hybracter_long/processing/pre_polish/BBR060204_cB0304_chromosome.fasta, results/BBR060204_cB0304/hybracter_long/processing/pre_polish/BBR060204_cB0304_ignore_list.txt
    output: results/BBR060204_cB0304/hybracter_long/processing/complete/dnaapler/BBR060204_cB0304/BBR060204_cB0304_reoriented.fasta, results/BBR060204_cB0304/hybracter_long/versions/BBR060204_cB0304/dnaapler.version
    log: results/BBR060204_cB0304/hybracter_long/stderr/dnaapler/BBR060204_cB0304.log
    jobid: 28
    benchmark: results/BBR060204_cB0304/hybracter_long/benchmarks/dnaapler/BBR060204_cB0304.txt
    reason: Missing output files: results/BBR060204_cB0304/hybracter_long/processing/complete/dnaapler/BBR060204_cB0304/BBR060204_cB0304_reoriented.fasta; Input files updated by another job: results/BBR060204_cB0304/hybracter_long/processing/pre_polish/BBR060204_cB0304_ignore_list.txt, results/BBR060204_cB0304/hybracter_long/processing/pre_polish/BBR060204_cB0304_chromosome.fasta; Code has changed since last execution; Params have changed since last execution
    wildcards: sample=BBR060204_cB0304
    threads: 8
    resources: tmpdir=/tmp, mem_mb=15259, mem_mib=15259, mem=16000MB, time=08:00:00


        dnaapler all -i results/BBR060204_cB0304/hybracter_long/processing/pre_polish/BBR060204_cB0304_chromosome.fasta -o results/BBR060204_cB0304/hybracter_long/processing/complete/dnaapler/BBR060204_cB0304 --ignore results/BBR060204_cB0304/hybracter_long/processing/pre_polish/BBR060204_cB0304_ignore_list.txt -p BBR060204_cB0304 -t 8 -a nearest --db dnaa,repa -f 2> results/BBR060204_cB0304/hybracter_long/stderr/dnaapler/BBR060204_cB0304.log
        dnaapler --version > results/BBR060204_cB0304/hybracter_long/versions/BBR060204_cB0304/dnaapler.version
```
(Notice the `--db dnaa,repa` parameter)

Which resulted in the following log from dnaapler:

```bash
2024-08-19 16:07:25.603 | INFO     | dnaapler.utils.util:begin_dnaapler:71 - You are using dnaapler version 0.8.0
2024-08-19 16:07:25.603 | INFO     | dnaapler.utils.util:begin_dnaapler:72 - Repository homepage is https://github.com/gbouras13/dnaapler
2024-08-19 16:07:25.603 | INFO     | dnaapler.utils.util:begin_dnaapler:73 - Written by George Bouras: george.bouras@adelaide.edu.au
2024-08-19 16:07:25.603 | INFO     | dnaapler.utils.util:begin_dnaapler:74 - Your input FASTA is results/BBR060204_cB0304/hybracter_long/processing/pre_polish/BBR060204_cB0304_chromosome.fasta
2024-08-19 16:07:25.604 | INFO     | dnaapler.utils.util:begin_dnaapler:75 - Your output directory  is results/BBR060204_cB0304/hybracter_long/processing/complete/dnaapler/BBR060204_cB0304
2024-08-19 16:07:25.604 | INFO     | dnaapler.utils.util:begin_dnaapler:76 - You have specified 8 threads to use with blastx
2024-08-19 16:07:25.604 | INFO     | dnaapler.utils.util:begin_dnaapler:77 - You have specified dnaA,repA gene(s) to reorient your sequence
2024-08-19 16:07:25.604 | INFO     | dnaapler.utils.util:check_blast_version:117 - Checking BLAST installation.
2024-08-19 16:07:25.694 | INFO     | dnaapler.utils.util:check_blast_version:137 - BLAST version found is v2.16.0.
2024-08-19 16:07:25.694 | INFO     | dnaapler.utils.util:check_blast_version:147 - BLAST version is ok.
2024-08-19 16:07:25.695 | INFO     | dnaapler.utils.util:check_pyrodigal_version:92 - Checking pyrodigal installation.
2024-08-19 16:07:25.695 | INFO     | dnaapler.utils.util:check_pyrodigal_version:103 - Pyrodigal version is v3.5.1
2024-08-19 16:07:25.695 | INFO     | dnaapler.utils.util:check_pyrodigal_version:104 - Pyrodigal version is ok.
2024-08-19 16:07:25.696 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --input results/BBR060204_cB0304/hybracter_long/processing/pre_polish/BBR060204_cB0304_chromosome.fasta
2024-08-19 16:07:25.696 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --output results/BBR060204_cB0304/hybracter_long/processing/complete/dnaapler/BBR060204_cB0304
2024-08-19 16:07:25.696 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --threads 8
2024-08-19 16:07:25.696 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --force True
2024-08-19 16:07:25.697 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --prefix BBR060204_cB0304
2024-08-19 16:07:25.697 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --evalue 1e-10
2024-08-19 16:07:25.697 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --autocomplete nearest
2024-08-19 16:07:25.697 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --seed_value 13
2024-08-19 16:07:25.698 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --ignore results/BBR060204_cB0304/hybracter_long/processing/pre_polish/BBR060204_cB0304_ignore_list.txt
2024-08-19 16:07:25.698 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --custom_db 
2024-08-19 16:07:25.698 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --db dnaa,repa
2024-08-19 16:07:25.698 | INFO     | dnaapler.utils.validation:validate_fasta_all:100 - Checking that the input file results/BBR060204_cB0304/hybracter_long/processing/pre_polish/BBR060204_cB0304_chromosome.fasta is in FASTA format and has at least 1 entry.
2024-08-19 16:07:25.724 | INFO     | dnaapler.utils.validation:validate_fasta_all:107 - results/BBR060204_cB0304/hybracter_long/processing/pre_polish/BBR060204_cB0304_chromosome.fasta file checked.
2024-08-19 16:07:25.734 | INFO     | dnaapler.utils.validation:validate_fasta_all:116 - results/BBR060204_cB0304/hybracter_long/processing/pre_polish/BBR060204_cB0304_chromosome.fasta has only one entry.
2024-08-19 16:07:25.739 | INFO     | dnaapler.utils.validation:check_evalue:187 - You have specified an evalue of 1e-10.
2024-08-19 16:07:25.739 | INFO     | dnaapler:all:930 - You have specified contigs to ignore in results/BBR060204_cB0304/hybracter_long/processing/pre_polish/BBR060204_cB0304_ignore_list.txt.
2024-08-19 16:07:25.740 | INFO     | dnaapler.utils.external_tools:run:53 - Started running blastx -db /home/rks/phd/ouroboros/.snakemake/conda/37009296e852720e3812295fad89acac_/lib/python3.12/site-packages/hybracter/workflow/conda/61381f58265b1d3b2f5a23d599a4e749_/lib/python3.12/site-packages/dnaapler/db/dnaA_repA_db -evalue 1e-10 -num_threads 8 -outfmt ' 6 qseqid qlen sseqid slen length qstart qend sstart send pident nident gaps mismatch evalue bitscore qseq sseq ' -out results/BBR060204_cB0304/hybracter_long/processing/complete/dnaapler/BBR060204_cB0304/BBR060204_cB0304_blast_output.txt -query results/BBR060204_cB0304/hybracter_long/processing/pre_polish/BBR060204_cB0304_chromosome.fasta ...
2024-08-19 16:07:42.844 | INFO     | dnaapler.utils.external_tools:run:55 - Done running blastx -db /home/rks/phd/ouroboros/.snakemake/conda/37009296e852720e3812295fad89acac_/lib/python3.12/site-packages/hybracter/workflow/conda/61381f58265b1d3b2f5a23d599a4e749_/lib/python3.12/site-packages/dnaapler/db/dnaA_repA_db -evalue 1e-10 -num_threads 8 -outfmt ' 6 qseqid qlen sseqid slen length qstart qend sstart send pident nident gaps mismatch evalue bitscore qseq sseq ' -out results/BBR060204_cB0304/hybracter_long/processing/complete/dnaapler/BBR060204_cB0304/BBR060204_cB0304_blast_output.txt -query results/BBR060204_cB0304/hybracter_long/processing/pre_polish/BBR060204_cB0304_chromosome.fasta
2024-08-19 16:07:42.845 | WARNING  | dnaapler:all:973 - results/BBR060204_cB0304/hybracter_long/processing/pre_polish/BBR060204_cB0304_ignore_list.txt contains no text. No contigs will be ignored
2024-08-19 16:07:42.860 | WARNING  | dnaapler.utils.processing:reorient_single_record_bulk:391 - The top blastx hit for the contig contig_1 did not begin with a valid start codon.
2024-08-19 16:07:42.860 | WARNING  | dnaapler.utils.processing:reorient_single_record_bulk:394 - Searching with pyrodigal for the CDS overlapping the most with the top blastx hit to reorient with.
2024-08-19 16:07:52.889 | INFO     | dnaapler.utils.util:end_dnaapler:160 - dnaapler has finished
2024-08-19 16:07:52.890 | INFO     | dnaapler.utils.util:end_dnaapler:161 - Elapsed time: 27.29 seconds
```

and the following output:

```bash
Contig	Gene_Reoriented	Start	Strand	Top_Hit	Top_Hit_Length	Covered_Length	Coverage	Identical_AAs	Identity_Percentage
contig_1	repA	1504958	forward	UniRef90_A0A0T8UAE4	460	236	51.3	81	34.32
```
Which incorrectly identifies a repA gene.

When I use the following command instead (which uses the same input as above, but removes the `--db` parameter):
```bash
dnaapler all -i results/BBR060204_cB0304/hybracter_long/processing/pre_polish/BBR060204_cB0304_chromosome.fasta -o BBR060204_cB0304_manual --ignore results/BBR060204_cB0304/hybracter_long/processing/pre_polish/BBR060204_cB0304_ignore_list.txt -p BBR060204_cB0304_manual -t 8 -a nearest -f
```

I get the following log:

```bash
2024-09-10 15:27:27.134 | INFO     | dnaapler.utils.util:begin_dnaapler:71 - You are using dnaapler version 0.8.0
2024-09-10 15:27:27.134 | INFO     | dnaapler.utils.util:begin_dnaapler:72 - Repository homepage is https://github.com/gbouras13/dnaapler
2024-09-10 15:27:27.135 | INFO     | dnaapler.utils.util:begin_dnaapler:73 - Written by George Bouras: george.bouras@adelaide.edu.au
2024-09-10 15:27:27.135 | INFO     | dnaapler.utils.util:begin_dnaapler:74 - Your input FASTA is results/BBR060204_cB0304/hybracter_long/processing/pre_polish/BBR060204_cB0304_chromosome.fasta
2024-09-10 15:27:27.135 | INFO     | dnaapler.utils.util:begin_dnaapler:75 - Your output directory  is BBR060204_cB0304_manual
2024-09-10 15:27:27.135 | INFO     | dnaapler.utils.util:begin_dnaapler:76 - You have specified 8 threads to use with blastx
2024-09-10 15:27:27.135 | INFO     | dnaapler.utils.util:begin_dnaapler:77 - You have specified all gene(s) to reorient your sequence
2024-09-10 15:27:27.135 | INFO     | dnaapler.utils.util:check_blast_version:117 - Checking BLAST installation.
2024-09-10 15:27:27.214 | INFO     | dnaapler.utils.util:check_blast_version:137 - BLAST version found is v2.16.0.
2024-09-10 15:27:27.214 | INFO     | dnaapler.utils.util:check_blast_version:147 - BLAST version is ok.
2024-09-10 15:27:27.215 | INFO     | dnaapler.utils.util:check_pyrodigal_version:92 - Checking pyrodigal installation.
2024-09-10 15:27:27.215 | INFO     | dnaapler.utils.util:check_pyrodigal_version:103 - Pyrodigal version is v3.5.1
2024-09-10 15:27:27.215 | INFO     | dnaapler.utils.util:check_pyrodigal_version:104 - Pyrodigal version is ok.
2024-09-10 15:27:27.215 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --input results/BBR060204_cB0304/hybracter_long/processing/pre_polish/BBR060204_cB0304_chromosome.fasta
2024-09-10 15:27:27.215 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --output BBR060204_cB0304_manual
2024-09-10 15:27:27.215 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --threads 8
2024-09-10 15:27:27.215 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --force True
2024-09-10 15:27:27.216 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --prefix BBR060204_cB0304_manual
2024-09-10 15:27:27.216 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --evalue 1e-10
2024-09-10 15:27:27.216 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --autocomplete nearest
2024-09-10 15:27:27.216 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --seed_value 13
2024-09-10 15:27:27.216 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --ignore results/BBR060204_cB0304/hybracter_long/processing/pre_polish/BBR060204_cB0304_ignore_list.txt
2024-09-10 15:27:27.216 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --custom_db 
2024-09-10 15:27:27.216 | INFO     | dnaapler.utils.util:begin_dnaapler:82 - Parameter: --db all
2024-09-10 15:27:27.216 | INFO     | dnaapler.utils.validation:validate_fasta_all:100 - Checking that the input file results/BBR060204_cB0304/hybracter_long/processing/pre_polish/BBR060204_cB0304_chromosome.fasta is in FASTA format and has at least 1 entry.
2024-09-10 15:27:27.229 | INFO     | dnaapler.utils.validation:validate_fasta_all:107 - results/BBR060204_cB0304/hybracter_long/processing/pre_polish/BBR060204_cB0304_chromosome.fasta file checked.
2024-09-10 15:27:27.236 | INFO     | dnaapler.utils.validation:validate_fasta_all:116 - results/BBR060204_cB0304/hybracter_long/processing/pre_polish/BBR060204_cB0304_chromosome.fasta has only one entry.
2024-09-10 15:27:27.241 | INFO     | dnaapler.utils.validation:check_evalue:187 - You have specified an evalue of 1e-10.
2024-09-10 15:27:27.241 | INFO     | dnaapler:all:930 - You have specified contigs to ignore in results/BBR060204_cB0304/hybracter_long/processing/pre_polish/BBR060204_cB0304_ignore_list.txt.
2024-09-10 15:27:27.242 | INFO     | dnaapler.utils.external_tools:run:53 - Started running blastx -db /home/rks/phd/ouroboros/.snakemake/conda/37009296e852720e3812295fad89acac_/lib/python3.12/site-packages/hybracter/workflow/conda/61381f58265b1d3b2f5a23d599a4e749_/lib/python3.12/site-packages/dnaapler/db/all_db -evalue 1e-10 -num_threads 8 -outfmt ' 6 qseqid qlen sseqid slen length qstart qend sstart send pident nident gaps mismatch evalue bitscore qseq sseq ' -out BBR060204_cB0304_manual/BBR060204_cB0304_manual_blast_output.txt -query results/BBR060204_cB0304/hybracter_long/processing/pre_polish/BBR060204_cB0304_chromosome.fasta ...
2024-09-10 15:29:15.087 | INFO     | dnaapler.utils.external_tools:run:55 - Done running blastx -db /home/rks/phd/ouroboros/.snakemake/conda/37009296e852720e3812295fad89acac_/lib/python3.12/site-packages/hybracter/workflow/conda/61381f58265b1d3b2f5a23d599a4e749_/lib/python3.12/site-packages/dnaapler/db/all_db -evalue 1e-10 -num_threads 8 -outfmt ' 6 qseqid qlen sseqid slen length qstart qend sstart send pident nident gaps mismatch evalue bitscore qseq sseq ' -out BBR060204_cB0304_manual/BBR060204_cB0304_manual_blast_output.txt -query results/BBR060204_cB0304/hybracter_long/processing/pre_polish/BBR060204_cB0304_chromosome.fasta
2024-09-10 15:29:15.088 | WARNING  | dnaapler:all:973 - results/BBR060204_cB0304/hybracter_long/processing/pre_polish/BBR060204_cB0304_ignore_list.txt contains no text. No contigs will be ignored
2024-09-10 15:29:15.121 | INFO     | dnaapler.utils.util:end_dnaapler:160 - dnaapler has finished
2024-09-10 15:29:15.121 | INFO     | dnaapler.utils.util:end_dnaapler:161 - Elapsed time: 107.99 seconds
```

With the following result:
```bash
Contig	Gene_Reoriented	Start	Strand	Top_Hit	Top_Hit_Length	Covered_Length	Coverage	Identical_AAs	Identity_Percentage
contig_1	cog1474	1409732	forward	cog1474_gi|341582567|ref|YP_004763059.1|	401	391	97.51	371	94.88
```

As you can see, when the `--db` parameter is ommited, the automatic selection is actually able to identify the correct cog1474 gene, as it has a much better Identity_Percentage.
Hi @richardstoeckl ,

Thanks so much for this. 

I think I will make a little modification to dnaapler so I can pass repA, dnaa and cog1474 (without terL to avoid reorienting bacteriophages). Then I will update hybracter.

Nonetheless, glad it worked and hopefully it is useful for you and other archaea researchers!

George



Fix patch:
diff --git a/.github/workflows/release.yml b/.github/workflows/release.yml
index eb350310..593dbf34 100644
--- a/.github/workflows/release.yml
+++ b/.github/workflows/release.yml
@@ -12,13 +12,13 @@ jobs:
 
     steps:
     - uses: actions/checkout@v2
-    - uses: conda-incubator/setup-miniconda@v2
+    - uses: conda-incubator/setup-miniconda@v3
       with:
         python-version: 3.11
         activate-environment: hybracterENV
         environment-file: build/environment.yaml
         auto-activate-base: false
-        miniforge-variant: Mambaforge
+        channel-priority: strict
         channels: conda-forge,bioconda,defaults
         auto-update-conda: true
     - name: Install dependencies
diff --git a/HISTORY.md b/HISTORY.md
index 1bd07345..9ba2d7fc 100644
--- a/HISTORY.md
+++ b/HISTORY.md
@@ -1,5 +1,37 @@
 # History
 
+
+### v0.9.0 Updates (18 September 2024)
+
+**`--auto` for automatic estimation of chromosome size**
+
+* Thanks to an [issue](https://github.com/gbouras13/hybracter/issues/90) and code from @[richardstoeckl](https://github.com/richardstoeckl), Hybracter can now estimate the estimated chromosome size for each sample by passing `--auto`. 
+* The implementation uses [kmc](https://github.com/refresh-bio/KMC). Specifically, Hybracter uses kmc to count the number of unique 21mers that appear at least 10 times in your long-read FASTQ file. This is because, for a given assembly of length L,  and a k-mer size of k, the total number of unique possible k-mers  will be given by ( L – k ) + 1, and if L >> k, then it suffices as an estimate of total assembly size
+* The estimated chromosome size used by Hybracter will actually be 80% of the number of 21-mers found at least 10 times, as it needs to account for plasmids
+* If you aren't sure whether you have enough data for assembly (i.e. coverage lower than 20x), be careful using `--auto`, because the actual assembly size will tend to be larger than the number of unique 21mers found at least 10 times. Therefore, the estimated chromosome size will almost certainly be an underestimate and may lead to Hybracter considering your assembly "complete" when in fact it isn't.
+
+* If you use `--auto`, you do not need to specify the chromosome length in the input. This means you don't need to `-c` with `long-single` or `hybrid-single` and in the input csv sample sheet, you do not need a column with chromosome length.
+  
+e.g. for `hybracter long` you only need 2 columns with sample name and long-read FASTQ file path:
+
+```bash
+s_aureus_sample1,sample1_long_read.fastq.gz
+p_aeruginosa_sample2,sample2_long_read.fastq.gz
+```
+
+and for `hybracter hybrid` you only need 4 columns with sample name, long-read FASTQ, and R1 and R2 short-read FASTQ file paths:
+
+```bash
+s_aureus_sample1,sample1_long_read.fastq.gz,sample1_SR_R1.fastq.gz,sample1_SR_R2.fastq.gz
+p_aeruginosa_sample2,sample2_long_read.fastq.gz,sample2_SR_R1.fastq.gz,sample2_SR_R2.fastq.gz
+```
+
+**Other changes**
+
+* Hybracter v0.9.0 will automatically support the reorientation of archaeal chromosomes (thanks @[richardstoeckl](https://github.com/richardstoeckl)) to begin with the cog1474 Orc1/cdc6 gene.
+* `--datadir` can now also accept 2 paths separated by a comma, if you have long reads and short reads in separate directories e.g. `--datadir "long_read_dir,short_read_dir"` (https://github.com/gbouras13/hybracter/issues/76).
+* `--min_depth` parameter added. Hybracter will error out if your QC'd long reads have a coverage lower than `min_depth` for a sample (https://github.com/gbouras13/hybracter/issues/89).
+
 ## v0.8.0 (30 August 2024)
 
 * Add `--datadir` that removes the need to add full paths in sample sheet (thanks @oschwengers)
diff --git a/README.md b/README.md
index 7789efc0..47efca9e 100644
--- a/README.md
+++ b/README.md
@@ -21,24 +21,22 @@
 - [Hybracter: Enabling Scalable, Automated, Complete and Accurate Bacterial Genome Assemblies](#hybracter-enabling-scalable-automated-complete-and-accurate-bacterial-genome-assemblies)
   - [Table of Contents](#table-of-contents)
   - [Quick Start](#quick-start)
-    - [Google Colab Notebooks](#google-colab-notebooks)
-    - [Mamba/Conda](#mambaconda)
+    - [Conda](#conda)
     - [Container](#container)
+    - [Google Colab Notebooks](#google-colab-notebooks)
   - [Documentation](#documentation)
   - [Manuscript](#manuscript)
   - [Description](#description)
   - [Pipeline](#pipeline)
   - [Benchmarking](#benchmarking)
-  - [v0.7.0 Updates (04 March 2024)](#v070-updates-04-march-2024)
-  - [v0.5.0 Updates (08 January 2024)](#v050-updates-08-january-2024)
-  - [v0.4.0 Updates (14 November 2023)](#v040-updates-14-november-2023)
-  - [v0.2.0 Updates 26 October 2023 - Medaka, Polishing and `--no_medaka`](#v020-updates-26-october-2023---medaka-polishing-and---no_medaka)
+  - [Recent Updates](#recent-updates)
+    - [v0.9.0 Updates (18 September 2024)](#v090-updates-18-september-2024)
   - [Why Would You Run Hybracter?](#why-would-you-run-hybracter)
   - [Other Options](#other-options)
       - [Trycycler](#trycycler)
       - [Dragonflye](#dragonflye)
   - [Installation](#installation)
-    - [Conda](#conda)
+    - [Conda](#conda-1)
     - [Pip](#pip)
     - [Source](#source)
   - [Main Commands](#main-commands)
@@ -56,6 +54,11 @@
     - [Main Output Files](#main-output-files)
   - [Snakemake Profiles](#snakemake-profiles)
   - [Advanced Configuration](#advanced-configuration)
+  - [Older Updates](#older-updates)
+    - [v0.7.0 Updates (04 March 2024)](#v070-updates-04-march-2024)
+    - [v0.5.0 Updates (08 January 2024)](#v050-updates-08-january-2024)
+    - [v0.4.0 Updates (14 November 2023)](#v040-updates-14-november-2023)
+    - [v0.2.0 Updates 26 October 2023 - Medaka, Polishing and `--no_medaka`](#v020-updates-26-october-2023---medaka-polishing-and---no_medaka)
   - [Version Log](#version-log)
   - [System](#system)
   - [Bugs and Suggestions](#bugs-and-suggestions)
@@ -64,13 +67,7 @@
 
 ## Quick Start
 
-### Google Colab Notebooks
-
-If you don't want to install `hybracter` locally, you can run it without any code using the colab notebook [https://colab.research.google.com/github/gbouras13/hybracter/blob/main/run_hybracter.ipynb](https://colab.research.google.com/github/gbouras13/hybracter/blob/main/run_hybracter.ipynb)
-
-This is only recommend if you have one or a few samples to assemble (it takes a while per sample due to the limited nature of Google Colab resources - probably an hour or two a sample). If you have more than this, a local install as described below is suggested.
-
-### Mamba/Conda
+### Conda
 
 `hybracter` is available to install with `pip` or `conda`.
 
@@ -79,13 +76,13 @@ You will need conda or mamba available so `hybracter` can install all the requir
 Therefore, it is recommended to install `hybracter` into a conda environment as follows.
 
 ```bash
-mamba create -n hybracterENV -c bioconda -c conda-forge  hybracter
+conda create -n hybracterENV -c bioconda -c conda-forge  hybracter
 conda activate hybracterENV
 hybracter --help
 hybracter install
 ```
 
-Mamba is **highly highly** recommend. Please see the [documentation](https://hybracter.readthedocs.io/en/latest/install/) for more details on how to install mamba.
+Miniforge is **highly highly** recommended. Please see the [documentation](https://hybracter.readthedocs.io/en/latest/install/) for more details on how to install miniforeg.
 
 When you run `hybracter` for the first time, all the required dependencies will be installed as required, so it will take longer than usual (usually a few minutes). Every time you run it afterwards, it will be a lot faster as the dependenices will be installed.
 
@@ -114,9 +111,15 @@ containerImage="$IMAGE_DIR/hybracter_0.8.0.sif"
 # example command with test fastqs
  singularity exec $containerImage    hybracter hybrid-single -l test_data/Fastqs/test_long_reads.fastq.gz \
  -1 test_data/Fastqs/test_short_reads_R1.fastq.gz  -2 test_data/Fastqs/test_short_reads_R2.fastq.gz \
- -o output_test_singularity -t 4 -c 50000
+ -o output_test_singularity -t 4 --auto
 ```
 
+### Google Colab Notebooks
+
+If you don't want to install `hybracter` locally, you can run it without any code using the colab notebook [https://colab.research.google.com/github/gbouras13/hybracter/blob/main/run_hybracter.ipynb](https://colab.research.google.com/github/gbouras13/hybracter/blob/main/run_hybracter.ipynb)
+
+This is only recommend if you have one or a few samples to assemble (it takes a while per sample due to the limited nature of Google Colab resources - probably an hour or two a sample). If you have more than this, a local install as described below is suggested.
+
 ## Documentation
 
 Documentation for `hybracter` is available [here](https://hybracter.readthedocs.io/en/latest/).
@@ -167,59 +170,38 @@ To summarise the conclusions:
 * If you want the fastest possible chromosome assemblies for applications like species ID or sequence typing that retain a high level of accuracy, Dragonflye is a good option.
 * Dragonflye should not be used if you care about recovering plasmids.
 
-## v0.7.0 Updates (04 March 2024)
-
-**Changes to short read polishing**
-
-* Logic added to run `polypolish` v0.6.0 with `--careful` and skip pypolca if the SR coverage estimate is below 5x (note: FASTA files for pypolca will be generated in the processing directory to play nice with Snakemake, but these will be identical to the polypolish output). 
-* For 5-25x coverage, `polypolish --careful` and `pypolca` with `--careful` will be run. 
-* For >25x coverage, `polypolish` default and `pypolca` with `--careful` will be run. 
-* A preprint justifying these changes will be available soon.
-
-**`--logic` changes**
-
-* By default, `--logic` defaults to `last` for `hybracter hybrid`, as there we have found that the polishing strategy implemented above never makes the assembly worse. We suggest never using `--logic best` with `hybracter hybrid`.
-
-**Changes for chromosome contigs and circularity**
-
-* If hybracter assembles a contig that is greater than the minimum chromosome length but not marked as circular by Flye, this will now be denoted as a chromosome, but not circular. The genome will be marked as complete also. 
-    * These will usually be assemblies with some issue (e.g. prophages, circularisation issues, heterogeneity) and probably require some more attention.
-    * For example, with the _Vibrio cholerae_ larger chromosome described [here](https://rrwick.github.io/2024/02/15/misassemblies.html), the genome will be marked as 'complete' but the contig will not be marked as 'circular' in the `hybracter` output.
-    * Such contigs will be polished and be in the final `_chromosome.fasta` output, but they will not be rotated by `dnaapler`. 
-    * These were previously being excluded, which was missing assemblies with structural heterogeneity (causing the chromosome not to completely circularise) or even bacteria with linear chromosomes like [_Borrelia_](https://www.nature.com/articles/37551). 
-
-**Adds `--depth_filter`** 
+## Recent Updates
 
-* This is passed to [Plassembler](https://github.com/gbouras13/plassembler) and will filter out all putative plasmid contigs that are lower than this depth fraction compared to the chromosome.
-* Defaults to 0.25 like Unicycler's implementation.
+### v0.9.0 Updates (18 September 2024)
 
-## v0.5.0 Updates (08 January 2024)
+**`--auto` for automatic estimation of chromosome size**
 
-Ryan Wick recently ran `hybracter long` on the latest Dorado v0.5.0 basecalled Nanopore reads (his [blog post](https://rrwick.github.io/2023/12/18/ont-only-accuracy-update.html)). You can read a write-up of the results [here](https://hybracter.readthedocs.io/en/latest/dorado_ryan_louise_0_5_0/). As a result, subsampling has been added to Hybracter. 
+* Thanks to an [issue](https://github.com/gbouras13/hybracter/issues/90) and code from @[richardstoeckl](https://github.com/richardstoeckl), Hybracter can now estimate the estimated chromosome size for each sample by passing `--auto`. 
+* The implementation uses [kmc](https://github.com/refresh-bio/KMC). Specifically, Hybracter uses kmc to count the number of unique 21mers that appear at least 10 times in your long-read FASTQ file. This is because, for a given assembly of length L,  and a k-mer size of k, the total number of unique possible k-mers  will be given by ( L – k ) + 1, and if L >> k, then it suffices as an estimate of total assembly size
+* The estimated chromosome size used by Hybracter will actually be 80% of the number of 21-mers found at least 10 times, as it needs to account for plasmids
+* If you aren't sure whether you have enough data for assembly (i.e. coverage lower than 20x), be careful using `--auto`, because the actual assembly size will tend to be larger than the number of unique 21mers found at least 10 times. Therefore, the estimated chromosome size will almost certainly be an underestimate and may lead to Hybracter considering your assembly "complete" when in fact it isn't.
 
-* Adds subsampling using `--subsample_depth` using Filtlong, based on some benchmarking of Dorado v0.5.0. Defaults to 100x of the estimated chromosome size `-c`.
-* Also adds stricter criteria for complete assemblies (aka ensures that identified chromosomes must be circularised according to Flye).
+* If you use `--auto`, you do not need to specify the chromosome length in the input. This means you don't need to `-c` with `long-single` or `hybrid-single` and in the input csv sample sheet, you do not need a column with chromosome length.
   
-## v0.4.0 Updates (14 November 2023)
+e.g. for `hybracter long` you only need 2 columns with sample name and long-read FASTQ file path:
 
-* Adds `--logic` parameter. You have 2 choices: `--logic best` (the default) or `--logic last`.
-* `--logic best` will run `hybracter` as normal and the best assembly (by ALE or pyrodigal mean length) will be selected as the final assembly.
-* `--logic last` will force hybracter to pick the last polished round as the final assembly even if it is not the best as per ALE/pyrodigal. So for `hybracter hybrid` this will default to the pypolca polished round, for hybracter long it will be Medaka round 2. You may wish to use this if you want all your isolates to be consistently assembled.
-* Adds reorientation of pre polished chromosome in case it is selected as the best assembly
-* Adds fixes to the chromosome comparisons - now it is much easier to interpret any changes between polishing rounds.
-
-
-## v0.2.0 Updates 26 October 2023 - Medaka, Polishing and `--no_medaka`
-
-Ryan Wick's [blogpost](https://rrwick.github.io/2023/10/24/ont-only-accuracy-update.html) on 24 October 2023 suggests that if you have new 5Hz SUP or Res (bacterial model specific) ONT reads, Medaka polishing often makes things worse! It also implies that Nanopore reads are almost good enough to assemble perfect bacterial genomes (at least with Trycycler) which is pretty awesome.
+```bash
+s_aureus_sample1,sample1_long_read.fastq.gz
+p_aeruginosa_sample2,sample2_long_read.fastq.gz
+```
 
-Combined with the difficulty and randomness in installing Medaka from Nanopore, I have therefore decided to add a `--no_medaka` flag into v0.2.0. 
+and for `hybracter hybrid` you only need 4 columns with sample name, long-read FASTQ, and R1 and R2 short-read FASTQ file paths:
 
-I have also set Medaka to be v1.8.0 and I do not intend to upgrade this going forward, as this is the most recent stable bioconda version that doesn't seem to cause too much grief. 
+```bash
+s_aureus_sample1,sample1_long_read.fastq.gz,sample1_SR_R1.fastq.gz,sample1_SR_R2.fastq.gz
+p_aeruginosa_sample2,sample2_long_read.fastq.gz,sample2_SR_R1.fastq.gz,sample2_SR_R2.fastq.gz
+```
 
-If you have trouble with Medaka installation, I'd therefore suggest please using `--no_medaka`.
+**Other changes**
 
-`hybracter` should still handle cases where Medaka makes assemblies worse. If Medaka makes your assembly appreciably worse, `hybracter` should choose the best most accurate assembly as the unpolished one in long mode. 
+* Hybracter v0.9.0 will automatically support the reorientation of archaeal chromosomes (thanks @[richardstoeckl](https://github.com/richardstoeckl)) to begin with the cog1474 Orc1/cdc6 gene.
+* `--datadir` can now also accept 2 paths separated by a comma, if you have long reads and short reads in separate directories e.g. `--datadir "long_read_dir,short_read_dir"` (https://github.com/gbouras13/hybracter/issues/76).
+* `--min_depth` parameter added. Hybracter will error out if your QC'd long reads have a coverage lower than `min_depth` for a sample (https://github.com/gbouras13/hybracter/issues/89).
 
 ## Why Would You Run Hybracter?
 
@@ -253,14 +235,14 @@ If you are looking for the best possible (manual) bacterial assembly for a singl
 
 ## Installation
 
-You will need conda and **highly recommended** mamba to run `hybracter`, because it is required for the installation of each compartmentalised environment (e.g. Flye will have its own environment). Please see the [documentation](https://hybracter.readthedocs.io/en/latest/install/) for more details on how to install mamba.
+You will need conda  (**highly recommended** through miniforge) to run `hybracter`, because it is required for the installation of each compartmentalised environment (e.g. Flye will have its own environment). Please see the [documentation](https://hybracter.readthedocs.io/en/latest/install/) for more details on how to install miniforge.
 
 ### Conda
 
 `hybracter` is available to install with `conda`. To install `hybracter` into a conda enviornment called `hybracterENV`:
 
 ```bash
-mamba create -n hybracterENV hybracter
+conda create -n hybracterENV hybracter
 conda activate hybracterENV
 hybracter --help
 hybracter install
@@ -270,11 +252,11 @@ hybracter install
 
 `hybracter` is available to install with `pip` . 
 
-You will also need conda or mamba available so `hybracter` can install all the required dependencies. Therefore, it is recommended to install `hybracter` into a conda environment as follows.
+You will also need conda available so `hybracter` can install all the required dependencies. Therefore, it is recommended to install `hybracter` into a conda environment as follows.
 
 
 ```bash
-mamba create -n hybracterENV pip
+conda create -n hybracterENV pip
 conda activate hybracterENV
 pip install hybracter
 hybracter --help
@@ -356,6 +338,17 @@ s_aureus_sample1,sample1_long_read.fastq.gz,2500000,sample1_SR_R1.fastq.gz,sampl
 p_aeruginosa_sample2,sample2_long_read.fastq.gz,5500000,sample2_SR_R1.fastq.gz,sample2_SR_R2.fastq.gz
 ```
 
+**Using `--auto`**
+
+* If you use `--auto`, you can remove the column with the chromosome length
+
+e.g.
+
+```bash
+s_aureus_sample1,sample1_long_read.fastq.gz,sample1_SR_R1.fastq.gz,sample1_SR_R2.fastq.gz
+p_aeruginosa_sample2,sample2_long_read.fastq.gz,sample2_SR_R1.fastq.gz,sample2_SR_R2.fastq.gz
+```
+
 #### `hybracter long`
 
 `hybracter long` also requires an input csv with no headers, but only 3 columns.
@@ -373,6 +366,15 @@ s_aureus_sample1,sample1_long_read.fastq.gz,2500000
 p_aeruginosa_sample2,sample2_long_read.fastq.gz,5500000
 ```
 
+**Using `--auto`**
+
+* If you use `--auto`, you can remove the column with the chromosome length
+
+```bash
+s_aureus_sample1,sample1_long_read.fastq.gz
+p_aeruginosa_sample2,sample2_long_read.fastq.gz
+```
+
 ## Usage
 
 #### `hybracter install`
@@ -502,6 +504,63 @@ hybracter hybrid --input <input.csv> --output <output_dir> --threads <threads> -
 
 Thanks to its Snakemake backend, you can modify resource requirements for each job contained within `hybracter` using the configuration file. A defauly can be created using the `hybracter config` command. This can make it even more efficient in server environment, as many jobs can be more efficiently parallelised than the default settings. For more information, please see the [documentation](https://hybracter.readthedocs.io/en/latest/configuration/) 
 
+## Older Updates 
+
+### v0.7.0 Updates (04 March 2024)
+
+**Changes to short read polishing**
+
+* Logic added to run `polypolish` v0.6.0 with `--careful` and skip pypolca if the SR coverage estimate is below 5x (note: FASTA files for pypolca will be generated in the processing directory to play nice with Snakemake, but these will be identical to the polypolish output). 
+* For 5-25x coverage, `polypolish --careful` and `pypolca` with `--careful` will be run. 
+* For >25x coverage, `polypolish` default and `pypolca` with `--careful` will be run. 
+* A preprint justifying these changes will be available soon.
+
+**`--logic` changes**
+
+* By default, `--logic` defaults to `last` for `hybracter hybrid`, as there we have found that the polishing strategy implemented above never makes the assembly worse. We suggest never using `--logic best` with `hybracter hybrid`.
+
+**Changes for chromosome contigs and circularity**
+
+* If hybracter assembles a contig that is greater than the minimum chromosome length but not marked as circular by Flye, this will now be denoted as a chromosome, but not circular. The genome will be marked as complete also. 
+    * These will usually be assemblies with some issue (e.g. prophages, circularisation issues, heterogeneity) and probably require some more attention.
+    * For example, with the _Vibrio cholerae_ larger chromosome described [here](https://rrwick.github.io/2024/02/15/misassemblies.html), the genome will be marked as 'complete' but the contig will not be marked as 'circular' in the `hybracter` output.
+    * Such contigs will be polished and be in the final `_chromosome.fasta` output, but they will not be rotated by `dnaapler`. 
+    * These were previously being excluded, which was missing assemblies with structural heterogeneity (causing the chromosome not to completely circularise) or even bacteria with linear chromosomes like [_Borrelia_](https://www.nature.com/articles/37551). 
+
+**Adds `--depth_filter`** 
+
+* This is passed to [Plassembler](https://github.com/gbouras13/plassembler) and will filter out all putative plasmid contigs that are lower than this depth fraction compared to the chromosome.
+* Defaults to 0.25 like Unicycler's implementation.
+
+### v0.5.0 Updates (08 January 2024)
+
+Ryan Wick recently ran `hybracter long` on the latest Dorado v0.5.0 basecalled Nanopore reads (his [blog post](https://rrwick.github.io/2023/12/18/ont-only-accuracy-update.html)). You can read a write-up of the results [here](https://hybracter.readthedocs.io/en/latest/dorado_ryan_louise_0_5_0/). As a result, subsampling has been added to Hybracter. 
+
+* Adds subsampling using `--subsample_depth` using Filtlong, based on some benchmarking of Dorado v0.5.0. Defaults to 100x of the estimated chromosome size `-c`.
+* Also adds stricter criteria for complete assemblies (aka ensures that identified chromosomes must be circularised according to Flye).
+  
+### v0.4.0 Updates (14 November 2023)
+
+* Adds `--logic` parameter. You have 2 choices: `--logic best` (the default) or `--logic last`.
+* `--logic best` will run `hybracter` as normal and the best assembly (by ALE or pyrodigal mean length) will be selected as the final assembly.
+* `--logic last` will force hybracter to pick the last polished round as the final assembly even if it is not the best as per ALE/pyrodigal. So for `hybracter hybrid` this will default to the pypolca polished round, for hybracter long it will be Medaka round 2. You may wish to use this if you want all your isolates to be consistently assembled.
+* Adds reorientation of pre polished chromosome in case it is selected as the best assembly
+* Adds fixes to the chromosome comparisons - now it is much easier to interpret any changes between polishing rounds.
+
+
+### v0.2.0 Updates 26 October 2023 - Medaka, Polishing and `--no_medaka`
+
+Ryan Wick's [blogpost](https://rrwick.github.io/2023/10/24/ont-only-accuracy-update.html) on 24 October 2023 suggests that if you have new 5Hz SUP or Res (bacterial model specific) ONT reads, Medaka polishing often makes things worse! It also implies that Nanopore reads are almost good enough to assemble perfect bacterial genomes (at least with Trycycler) which is pretty awesome.
+
+Combined with the difficulty and randomness in installing Medaka from Nanopore, I have therefore decided to add a `--no_medaka` flag into v0.2.0. 
+
+I have also set Medaka to be v1.8.0 and I do not intend to upgrade this going forward, as this is the most recent stable bioconda version that doesn't seem to cause too much grief. 
+
+If you have trouble with Medaka installation, I'd therefore suggest please using `--no_medaka`.
+
+`hybracter` should still handle cases where Medaka makes assemblies worse. If Medaka makes your assembly appreciably worse, `hybracter` should choose the best most accurate assembly as the unpolished one in long mode. 
+
+
 ## Version Log
 
 A brief description of what is new in each update of `hybracter` can be found in the HISTORY.md file.
@@ -570,3 +629,5 @@ Pypolca:
 Snakemake:
 * Mölder F, Jablonski KP, Letcher B et al. Sustainable data analysis with Snakemake [version 1; peer review: 1 approved, 1 approved with reservations]. F1000Research 2021, 10:33 (https://doi.org/10.12688/f1000research.29032.1).
 
+KMC:
+* Marek Kokot, Maciej Długosz, Sebastian Deorowicz, KMC 3: counting and manipulating k-mer statistics, Bioinformatics, Volume 33, Issue 17, 01 September 2017, Pages 2759–2761, (https://doi.org/10.1093/bioinformatics/btx304).
\ No newline at end of file
diff --git a/docs/citation.md b/docs/citation.md
index 8b5c1ea9..587870e6 100644
--- a/docs/citation.md
+++ b/docs/citation.md
@@ -54,3 +54,5 @@ Pypolca:
 Snakemake:
 * Mölder F, Jablonski KP, Letcher B et al. Sustainable data analysis with Snakemake [version 1; peer review: 1 approved, 1 approved with reservations]. F1000Research 2021, 10:33 (https://doi.org/10.12688/f1000research.29032.1).
 
+KMC:
+* Marek Kokot, Maciej Długosz, Sebastian Deorowicz, KMC 3: counting and manipulating k-mer statistics, Bioinformatics, Volume 33, Issue 17, 01 September 2017, Pages 2759–2761, (https://doi.org/10.1093/bioinformatics/btx304).
diff --git a/docs/index.md b/docs/index.md
index 35abc1e5..d41c0eaf 100644
--- a/docs/index.md
+++ b/docs/index.md
@@ -13,7 +13,7 @@ It scales massively using the embarassingly parallel power of HPC and Snakemake
 
 ![Image](hybracter.png)
 
-- A. Reads are quality controlled and subsampled with [Filtlong](https://github.com/rrwick/Filtlong), [Porechop](https://github.com/rrwick/Porechop), [fastp](https://github.com/OpenGene/fastp), [Seqkit](https://github.com/shenwei356/seqkit) and optionally contaminant removal using modules from [trimnami](https://github.com/beardymcjohnface/Trimnami).
+- A. Reads are quality controlled and subsampled with [Filtlong](https://github.com/rrwick/Filtlong), [Porechop](https://github.com/rrwick/Porechop), [fastp](https://github.com/OpenGene/fastp), [Seqkit](https://github.com/shenwei356/seqkit) and optionally contaminant removal using modules from [trimnami](https://github.com/beardymcjohnface/Trimnami). Chromosome length is optionally estimated with [kmc](https://github.com/refresh-bio/KMC).
 - B. Long-read assembly is conducted with [Flye](https://github.com/fenderglass/Flye). Each sample is clssified if the chromosome(s) were assembled (marked as 'complete') or not (marked as 'incomplete') based on the given minimum chromosome length.
 - C. For complete isolates, plasmid recovery with [Plassembler](https://github.com/gbouras13/plassembler).
 - D. For all isolates, long read polishing with [Medaka](https://github.com/nanoporetech/medaka).
diff --git a/docs/input.md b/docs/input.md
index 0e8c6294..6de46edd 100644
--- a/docs/input.md
+++ b/docs/input.md
@@ -75,6 +75,21 @@ p_aeruginosa_sample2,sample2_long_read.fastq.gz,5500000,sample2_SR_R1.fastq.gz,s
 | -------------------- | -------------------------- | ------- | ---------------------- | ----------------------- |
 | p_aeruginosa_sample2 | sample2_long_read.fastq.gz | 5500000 | sample2_SR_R1.fastq.gz |  sample2_SR_R2.fastq.gz |
 
+
+* If you use `--auto`, you can remove the column with the chromosome length
+
+e.g.
+
+```bash
+s_aureus_sample1,sample1_long_read.fastq.gz,sample1_SR_R1.fastq.gz,sample1_SR_R2.fastq.gz
+p_aeruginosa_sample2,sample2_long_read.fastq.gz,sample2_SR_R1.fastq.gz,sample2_SR_R2.fastq.gz
+```
+
+| s_aureus_sample1     | sample1_long_read.fastq.gz | sample1_SR_R1.fastq.gz |  sample1_SR_R2.fastq.gz |
+| -------------------- | -------------------------- | ---------------------- | ----------------------- |
+| p_aeruginosa_sample2 | sample2_long_read.fastq.gz | sample2_SR_R1.fastq.gz |  sample2_SR_R2.fastq.gz |
+
+
 #### `hybracter long`
 
 `hybracter long` also requires an input csv with no headers, but only 3 columns.
@@ -96,6 +111,16 @@ p_aeruginosa_sample2,sample2_long_read.fastq.gz,5500000
 | -------------------- | -------------------------- | ------- |
 | p_aeruginosa_sample2 | sample2_long_read.fastq.gz | 5500000 |
 
+* If you use `--auto`, you can remove the column with the chromosome length
+
+```bash
+s_aureus_sample1,sample1_long_read.fastq.gz
+p_aeruginosa_sample2,sample2_long_read.fastq.gz
+```
+
+| s_aureus_sample1     | sample1_long_read.fastq.gz |
+| -------------------- | -------------------------- |
+| p_aeruginosa_sample2 | sample2_long_read.fastq.gz |
 
 ## Pacbio Reads
 
diff --git a/docs/install.md b/docs/install.md
index 5af818c4..abb74a39 100644
--- a/docs/install.md
+++ b/docs/install.md
@@ -8,14 +8,14 @@ This is only recommend if you have one or a few samples to assemble (it takes a
 
 # Local Installation
 
-You will need conda and **highly recommended** mamba to run `hybracter`, because it is required for the installation of each compartmentalised environment (e.g. Flye will have its own environment). See the end of this page for steps on how to install mamba.
+You will need conda, ideal via **highly recommended** [miniforge](https://github.com/conda-forge/miniforge) to run `hybracter`, because it is required for the installation of each compartmentalised environment (e.g. Flye will have its own environment). See the end of this page for steps on how to install [miniforge](https://github.com/conda-forge/miniforge).
 
 ## Conda
 
 `hybracter` is available to install with `conda`. To install `hybracter` into a conda enviornment called `hybracterENV`:
 
 ```bash
-mamba create -n hybracterENV  -c bioconda -c conda-forge  hybracter
+conda create -n hybracterENV  -c bioconda -c conda-forge  hybracter
 conda activate hybracterENV
 hybracter --help
 hybracter install
@@ -25,19 +25,17 @@ hybracter install
 
 `hybracter` is available to install with `pip` . 
 
-You will also need conda or mamba available so `hybracter` can install all the required dependencies. Therefore, it is recommended to install `hybracter` into a conda environment as follows.
+You will also need conda available so `hybracter` can install all the required dependencies. Therefore, it is recommended to install `hybracter` into a conda environment as follows.
 
 
 ```bash
-mamba create -n hybracterENV pip
+conda create -n hybracterENV pip
 conda activate hybracterENV
 pip install hybracter
 hybracter --help
 hybracter install
 ```
 
-Mamba is **highly highly** recommend. Please see the [documentation](https://hybracter.readthedocs.io/en/latest/install/) for more details on how to install mamba.
-
 ## Source
 
 Alternatively, the development version of `hybracter` (which may include new, untested features) can be installed manually via github. 
@@ -100,7 +98,7 @@ The containerised environments will automatically be created as required when yo
 
 For users who will be running `hybracter` offline e.g. on a cluster, this means you should probably run a small test datasets like `hybracter test-hybrid` and `hybracter test-long` where internet is available online (such as on the head node) to ensure all required environments are installed. I am one of these people :)
 
-Finally, but default `hybracter` will use mamba to install your environments. If for some reason you must use conda not mamba, use `--conda-frontend conda` with your hybracter command to force `hybracter` to use conda. 
+Finally, by default `hybracter` will use `conda` to install your environments (from v0.9.0), now that by default it uses the mamba solver. If you must use `mamba`, use `--conda-frontend mamba` with your hybracter command to force `hybracter` to use mamba explicitly. 
 
 ## Errors with Installing Dependencies & Environments
 
@@ -148,9 +146,9 @@ hybracter test-long --threads 1
 
 # Beginner Conda and Mamba Installation
 
-If you are new to using the command-line, please install mamba or conda using either option detailed below.
+If you are new to using the command-line, please install conda using either option detailed below.
 
-I would prefer and recommend the miniforge option as it comes with mamba!
+I would prefer and recommend the miniforge option!
 
 ##  Install miniforge
 
@@ -180,7 +178,7 @@ conda config --add channels conda-forge
 ##### Finally, I would recommend installing hybracter into a fresh environment. For example to create an environment called hybracterENV with hybracter installed:
 
 ```bash
-mamba create -n hybracterENV hybracter
+conda create -n hybracterENV hybracter
 conda activate hybracterENV
 hybracter -h
 ```
@@ -210,16 +208,12 @@ conda config --add channels bioconda
 conda config --add channels conda-forge
 ```
 
-##### After this, conda should be installed (you may need to restart your terminal). It is recommended that mamba is also installed, as it will solve the enviroment quicker than conda:
+##### I would recommend installing hybracter into a fresh environment. For example to create an environment called hybracterENV with hybracter installed:
 
 ```bash
-conda install mamba
-```
-
-##### Finally, I would recommend installing hybracter into a fresh environment. For example to create an environment called hybracterENV with hybracter installed:
-
-```bash
-mamba create -n hybracterENV hybracter
+conda create -n hybracterENV hybracter
 conda activate hybracterENV
 hybracter -h
 ```
+
+**There is no need to install mamba anymore, as conda's default solver is mamba!**
\ No newline at end of file
diff --git a/docs/output.md b/docs/output.md
index e7dfbc9e..5c8a3c84 100644
--- a/docs/output.md
+++ b/docs/output.md
@@ -93,7 +93,7 @@ The `processing` directory will contain a number of intermediate directories who
 
 ##### 2. `qc` directory
 
- This directory will contain the filtered, trimmed and contaminant removed FASTQ reads (where applicable) for each sample
+ This directory will contain the filtered, trimmed and contaminant removed FASTQ reads (where applicable) for each sample, along with the kmc estimated chromosome size (which will be used if you specify `--auto`)
 
 ##### 3. `plassembler` directory
 
diff --git a/docs/run.md b/docs/run.md
index 65662b43..f35d502a 100644
--- a/docs/run.md
+++ b/docs/run.md
@@ -37,12 +37,13 @@ hybracter hybrid -i <input.csv> -o <output_dir> -t <threads>  [other arguments]
 * You can turn off Medaka polishing using `--no_medaka` - recommended for Q20+ modern Nanopore reads
 * You can turn off pypolca polishing using `--no_pypolca` - I wouldn't though!
 * You can change the `--depth_filter` from 0.25x chromosome coverage. This will filter out all Plassembler contigs below this depth.
-* By default, `hybracter hybrid` takes the last polishing round as the final assembly  (`--logic last`). We would not recommend changing this to `--logic best`, as picking the best polishing round according to ALE with  `--logic best` is not guaranteed to give the most accurate assembly (See our [preprint](https://www.biorxiv.org/content/10.1101/2024.03.07.584013v1)).
+* By default, `hybracter hybrid` takes the last polishing round as the final assembly  (`--logic last`). We would not recommend changing this to `--logic best`, as picking the best polishing round according to ALE with  `--logic best` is not guaranteed to give the most accurate assembly (See our [paper](https://doi.org/10.1099/mgen.0.001244)).
+* You can estimate the chromosome size with kmc by using `--auto`
+* You can set a minimum long-read depth with `--min_depth`. Hybracter will error out if your estimated long-reads coverage is lower than this.
+* If you are running hybracter on a Mac, please use `--mac` (or find a Linux machine). This will make sure Medaka v1.8.0 is installed, as newer versions don't work on Macs.
 
 ```bash
-
-hybracter version 0.7.0
-
+hybracter version 0.9.0
 
 
  _           _                    _            
@@ -59,13 +60,19 @@ Usage: hybracter hybrid [OPTIONS] [SNAKE_ARGS]...
 
 Options:
   -i, --input TEXT                Input csv  [required]
+  --datadir TEXT                  Directory/ies where FASTQs are. Can specify
+                                  1 directory (long and short FASTQs in the
+                                  same directory) or 2 (long and short FASTQs
+                                  in separate directories). If you specify 2,
+                                  they must be separated by a comma e.g.
+                                  dirlong,dirshort. Will be added to the
+                                  filename if an input csv if provided.
   --no_pypolca                    Do not use pypolca to polish assemblies with
                                   short reads
   --logic [best|last]             Hybracter logic to select best assembly. Use
+                                  --last to pick the last polishing round. Use
                                   --best to pick best assembly based on ALE
-                                  (hybrid) or pyrodigal mean length (long).
-                                  Use --last to pick the last polishing round
-                                  regardless.  [default: last]
+                                  (hybrid).   [default: last]
   -o, --output PATH               Output directory  [default: hybracter_out]
   --configfile TEXT               Custom config file [default: config.yaml]
   -t, --threads INTEGER           Number of threads to use  [default: 1]
@@ -79,9 +86,15 @@ Options:
   --subsample_depth INTEGER       subsampled long read depth to subsample with
                                   Filtlong. By default is 100x.  [default:
                                   100]
-  --medakaModel [r1041_e82_400bps_hac_v4.2.0|r1041_e82_400bps_sup_v4.2.0|r941_sup_plant_g610|r941_min_fast_g507|r941_prom_fast_g507|r941_min_fast_g303|r941_min_high_g303|r941_min_high_g330|r941_prom_fast_g303|r941_prom_high_g303|r941_prom_high_g330|r941_min_high_g344|r941_min_high_g351|r941_min_high_g360|r941_prom_high_g344|r941_prom_high_g360|r941_prom_high_g4011|r10_min_high_g303|r10_min_high_g340|r103_min_high_g345|r103_min_high_g360|r103_prom_high_g360|r103_fast_g507|r103_hac_g507|r103_sup_g507|r104_e81_fast_g5015|r104_e81_sup_g5015|r104_e81_hac_g5015|r104_e81_sup_g610|r1041_e82_400bps_hac_g615|r1041_e82_400bps_fast_g615|r1041_e82_400bps_fast_g632|r1041_e82_260bps_fast_g632|r1041_e82_400bps_hac_g632|r1041_e82_400bps_sup_g615|r1041_e82_260bps_hac_g632|r1041_e82_260bps_sup_g632|r1041_e82_400bps_hac_v4.0.0|r1041_e82_400bps_sup_v4.0.0|r1041_e82_260bps_hac_v4.0.0|r1041_e82_260bps_sup_v4.0.0|r1041_e82_260bps_hac_v4.1.0|r1041_e82_260bps_sup_v4.1.0|r1041_e82_400bps_hac_v4.1.0|r1041_e82_400bps_sup_v4.1.0|r941_min_high_g340_rle|r941_min_hac_g507|r941_min_sup_g507|r941_prom_hac_g507|r941_prom_sup_g507|r941_e81_fast_g514|r941_e81_hac_g514|r941_e81_sup_g514]
+  --min_depth INTEGER             minimum long read depth to continue the run.
+                                  By default is 0x. Hybracter will error and
+                                  exit if a sample has less than
+                                  min_depth*chromosome_size bases of long-
+                                  reads left AFTER filtlong and porechop-ABI
+                                  steps are run.  [default: 0]
+  --medakaModel [r1041_e82_400bps_sup_v5.0.0|r1041_e82_400bps_hac_v5.0.0|r1041_e82_400bps_hac_v4.3.0|r1041_e82_400bps_sup_v4.3.0|r1041_e82_400bps_hac_v4.2.0|r1041_e82_400bps_sup_v4.2.0|r941_sup_plant_g610|r941_min_fast_g507|r941_prom_fast_g507|r941_min_fast_g303|r941_min_high_g303|r941_min_high_g330|r941_prom_fast_g303|r941_prom_high_g303|r941_prom_high_g330|r941_min_high_g344|r941_min_high_g351|r941_min_high_g360|r941_prom_high_g344|r941_prom_high_g360|r941_prom_high_g4011|r10_min_high_g303|r10_min_high_g340|r103_min_high_g345|r103_min_high_g360|r103_prom_high_g360|r103_fast_g507|r103_hac_g507|r103_sup_g507|r104_e81_fast_g5015|r104_e81_sup_g5015|r104_e81_hac_g5015|r104_e81_sup_g610|r1041_e82_400bps_hac_g615|r1041_e82_400bps_fast_g615|r1041_e82_400bps_fast_g632|r1041_e82_260bps_fast_g632|r1041_e82_400bps_hac_g632|r1041_e82_400bps_sup_g615|r1041_e82_260bps_hac_g632|r1041_e82_260bps_sup_g632|r1041_e82_400bps_hac_v4.0.0|r1041_e82_400bps_sup_v4.0.0|r1041_e82_260bps_hac_v4.0.0|r1041_e82_260bps_sup_v4.0.0|r1041_e82_260bps_hac_v4.1.0|r1041_e82_260bps_sup_v4.1.0|r1041_e82_400bps_hac_v4.1.0|r1041_e82_400bps_sup_v4.1.0|r941_min_high_g340_rle|r941_min_hac_g507|r941_min_sup_g507|r941_prom_hac_g507|r941_prom_sup_g507|r941_e81_fast_g514|r941_e81_hac_g514|r941_e81_sup_g514]
                                   Medaka Model.  [default:
-                                  r1041_e82_400bps_sup_v4.2.0]
+                                  r1041_e82_400bps_sup_v5.0.0]
   --flyeModel [--nano-hq|--nano-corr|--nano-raw|--pacbio-raw|--pacbio-corr|--pacbio-hifi]
                                   Flye Assembly Parameter  [default: --nano-
                                   hq]
@@ -93,11 +106,16 @@ Options:
                                   be used as a database with dnaapler custom.
   --no_medaka                     Do not polish the long read assembly with
                                   Medaka.
+  --auto                          Automatically estimate the chromosome size
+                                  using KMC.
   --depth_filter FLOAT            Depth filter to pass to Plassembler. Filters
                                   out all putative plasmid contigs below this
                                   fraction of the chromosome read depth (needs
                                   to be below in both long and short read sets
                                   for hybrid).
+  --mac                           If you are running Hybracter on Mac -
+                                  installs v1.8.0 of Medaka as higher versions
+                                  break.
   --use-conda / --no-use-conda    Use conda for Snakemake rules  [default:
                                   use-conda]
   --conda-prefix PATH             Custom conda env directory
@@ -106,12 +124,11 @@ Options:
                                   --nolock, --show-failed-logs, --conda-
                                   frontend mamba]
   -h, --help                      Show this message and exit.
-
 ```
 
 ## `hybracter hybrid-single`
 
-You can also run a single isolate using the same input arguments as Unicycler by specifying `hybracter hybrid-single`. Instead of specifying a CSV with `--input`, use `-l` to specify the long read FASTQ file, `-1` to specify the short read R1 file, `-2` to specify the short read R2 file, `-s` to specify the sample name, `-c` to specify the chromosome size.
+You can also run a single isolate using the same input arguments as Unicycler by specifying `hybracter hybrid-single`. Instead of specifying a CSV with `--input`, use `-l` to specify the long read FASTQ file, `-1` to specify the short read R1 file, `-2` to specify the short read R2 file, `-s` to specify the sample name, `-c` to specify the chromosome size (`-c` can be omitted with `--auto`).
 
 ```bash
 hybracter hybrid-single -l <longread FASTQ> -1 <R1 short reads FASTQ> -2 <R2 short reads FASTQ> -s <sample name> -c <chromosome size> -o <output_dir> -t <threads>  [other arguments]
@@ -154,9 +171,15 @@ Options:
   --subsample_depth INTEGER       subsampled long read depth to subsample with
                                   Filtlong. By default is 100x.  [default:
                                   100]
-  --medakaModel [r1041_e82_400bps_hac_v4.2.0|r1041_e82_400bps_sup_v4.2.0|r941_sup_plant_g610|r941_min_fast_g507|r941_prom_fast_g507|r941_min_fast_g303|r941_min_high_g303|r941_min_high_g330|r941_prom_fast_g303|r941_prom_high_g303|r941_prom_high_g330|r941_min_high_g344|r941_min_high_g351|r941_min_high_g360|r941_prom_high_g344|r941_prom_high_g360|r941_prom_high_g4011|r10_min_high_g303|r10_min_high_g340|r103_min_high_g345|r103_min_high_g360|r103_prom_high_g360|r103_fast_g507|r103_hac_g507|r103_sup_g507|r104_e81_fast_g5015|r104_e81_sup_g5015|r104_e81_hac_g5015|r104_e81_sup_g610|r1041_e82_400bps_hac_g615|r1041_e82_400bps_fast_g615|r1041_e82_400bps_fast_g632|r1041_e82_260bps_fast_g632|r1041_e82_400bps_hac_g632|r1041_e82_400bps_sup_g615|r1041_e82_260bps_hac_g632|r1041_e82_260bps_sup_g632|r1041_e82_400bps_hac_v4.0.0|r1041_e82_400bps_sup_v4.0.0|r1041_e82_260bps_hac_v4.0.0|r1041_e82_260bps_sup_v4.0.0|r1041_e82_260bps_hac_v4.1.0|r1041_e82_260bps_sup_v4.1.0|r1041_e82_400bps_hac_v4.1.0|r1041_e82_400bps_sup_v4.1.0|r941_min_high_g340_rle|r941_min_hac_g507|r941_min_sup_g507|r941_prom_hac_g507|r941_prom_sup_g507|r941_e81_fast_g514|r941_e81_hac_g514|r941_e81_sup_g514]
+  --min_depth INTEGER             minimum long read depth to continue the run.
+                                  By default is 0x. Hybracter will error and
+                                  exit if a sample has less than
+                                  min_depth*chromosome_size bases of long-
+                                  reads left AFTER filtlong and porechop-ABI
+                                  steps are run.  [default: 0]
+  --medakaModel [r1041_e82_400bps_sup_v5.0.0|r1041_e82_400bps_hac_v5.0.0|r1041_e82_400bps_hac_v4.3.0|r1041_e82_400bps_sup_v4.3.0|r1041_e82_400bps_hac_v4.2.0|r1041_e82_400bps_sup_v4.2.0|r941_sup_plant_g610|r941_min_fast_g507|r941_prom_fast_g507|r941_min_fast_g303|r941_min_high_g303|r941_min_high_g330|r941_prom_fast_g303|r941_prom_high_g303|r941_prom_high_g330|r941_min_high_g344|r941_min_high_g351|r941_min_high_g360|r941_prom_high_g344|r941_prom_high_g360|r941_prom_high_g4011|r10_min_high_g303|r10_min_high_g340|r103_min_high_g345|r103_min_high_g360|r103_prom_high_g360|r103_fast_g507|r103_hac_g507|r103_sup_g507|r104_e81_fast_g5015|r104_e81_sup_g5015|r104_e81_hac_g5015|r104_e81_sup_g610|r1041_e82_400bps_hac_g615|r1041_e82_400bps_fast_g615|r1041_e82_400bps_fast_g632|r1041_e82_260bps_fast_g632|r1041_e82_400bps_hac_g632|r1041_e82_400bps_sup_g615|r1041_e82_260bps_hac_g632|r1041_e82_260bps_sup_g632|r1041_e82_400bps_hac_v4.0.0|r1041_e82_400bps_sup_v4.0.0|r1041_e82_260bps_hac_v4.0.0|r1041_e82_260bps_sup_v4.0.0|r1041_e82_260bps_hac_v4.1.0|r1041_e82_260bps_sup_v4.1.0|r1041_e82_400bps_hac_v4.1.0|r1041_e82_400bps_sup_v4.1.0|r941_min_high_g340_rle|r941_min_hac_g507|r941_min_sup_g507|r941_prom_hac_g507|r941_prom_sup_g507|r941_e81_fast_g514|r941_e81_hac_g514|r941_e81_sup_g514]
                                   Medaka Model.  [default:
-                                  r1041_e82_400bps_sup_v4.2.0]
+                                  r1041_e82_400bps_sup_v5.0.0]
   --flyeModel [--nano-hq|--nano-corr|--nano-raw|--pacbio-raw|--pacbio-corr|--pacbio-hifi]
                                   Flye Assembly Parameter  [default: --nano-
                                   hq]
@@ -168,11 +191,16 @@ Options:
                                   be used as a database with dnaapler custom.
   --no_medaka                     Do not polish the long read assembly with
                                   Medaka.
+  --auto                          Automatically estimate the chromosome size
+                                  using KMC.
   --depth_filter FLOAT            Depth filter to pass to Plassembler. Filters
                                   out all putative plasmid contigs below this
                                   fraction of the chromosome read depth (needs
                                   to be below in both long and short read sets
                                   for hybrid).
+  --mac                           If you are running Hybracter on Mac -
+                                  installs v1.8.0 of Medaka as higher versions
+                                  break.
   --use-conda / --no-use-conda    Use conda for Snakemake rules  [default:
                                   use-conda]
   --conda-prefix PATH             Custom conda env directory
@@ -205,16 +233,18 @@ hybracter long -i <input.csv> -o <output_dir> -t <threads> [other arguments]
 * You can turn off Medaka polishing using `--no_medaka` - recommended for Q20+ modern Nanopore and PacBio reads
 * You can change the `--depth_filter` from 0.25x chromosome coverage. This will filter out all Plassembler contigs below this depth.
 * You can force `hybracter long` to pick the last polishing round (not the best according to pyrodigal mean CDS length) with `--logic last`. `hybracter long` defaults to picking the best i.e. `--logic best`.
-
-
+* You can estimate the chromosome size with kmc by using `--auto`
+* You can set a minimum long-read depth with `--min_depth`. Hybracter will error out if your estimated long-reads coverage is lower than this.
+* If you are running hybracter on a Mac, please use `--mac` (or find a Linux machine). This will make sure Medaka v1.8.0 is installed, as newer versions don't work on Macs.
 
 ```bash
 Usage: hybracter long [OPTIONS] [SNAKE_ARGS]...
 
   Run hybracter with only long reads
-
 Options:
   -i, --input TEXT                Input csv  [required]
+  --datadir TEXT                  Directory where FASTQs are. Will be added to
+                                  the filename if in input csv if provided
   -o, --output PATH               Output directory  [default: hybracter_out]
   --configfile TEXT               Custom config file [default: config.yaml]
   -t, --threads INTEGER           Number of threads to use  [default: 1]
@@ -228,9 +258,15 @@ Options:
   --subsample_depth INTEGER       subsampled long read depth to subsample with
                                   Filtlong. By default is 100x.  [default:
                                   100]
-  --medakaModel [r1041_e82_400bps_hac_v4.2.0|r1041_e82_400bps_sup_v4.2.0|r941_sup_plant_g610|r941_min_fast_g507|r941_prom_fast_g507|r941_min_fast_g303|r941_min_high_g303|r941_min_high_g330|r941_prom_fast_g303|r941_prom_high_g303|r941_prom_high_g330|r941_min_high_g344|r941_min_high_g351|r941_min_high_g360|r941_prom_high_g344|r941_prom_high_g360|r941_prom_high_g4011|r10_min_high_g303|r10_min_high_g340|r103_min_high_g345|r103_min_high_g360|r103_prom_high_g360|r103_fast_g507|r103_hac_g507|r103_sup_g507|r104_e81_fast_g5015|r104_e81_sup_g5015|r104_e81_hac_g5015|r104_e81_sup_g610|r1041_e82_400bps_hac_g615|r1041_e82_400bps_fast_g615|r1041_e82_400bps_fast_g632|r1041_e82_260bps_fast_g632|r1041_e82_400bps_hac_g632|r1041_e82_400bps_sup_g615|r1041_e82_260bps_hac_g632|r1041_e82_260bps_sup_g632|r1041_e82_400bps_hac_v4.0.0|r1041_e82_400bps_sup_v4.0.0|r1041_e82_260bps_hac_v4.0.0|r1041_e82_260bps_sup_v4.0.0|r1041_e82_260bps_hac_v4.1.0|r1041_e82_260bps_sup_v4.1.0|r1041_e82_400bps_hac_v4.1.0|r1041_e82_400bps_sup_v4.1.0|r941_min_high_g340_rle|r941_min_hac_g507|r941_min_sup_g507|r941_prom_hac_g507|r941_prom_sup_g507|r941_e81_fast_g514|r941_e81_hac_g514|r941_e81_sup_g514]
+  --min_depth INTEGER             minimum long read depth to continue the run.
+                                  By default is 0x. Hybracter will error and
+                                  exit if a sample has less than
+                                  min_depth*chromosome_size bases of long-
+                                  reads left AFTER filtlong and porechop-ABI
+                                  steps are run.  [default: 0]
+  --medakaModel [r1041_e82_400bps_sup_v5.0.0|r1041_e82_400bps_hac_v5.0.0|r1041_e82_400bps_hac_v4.3.0|r1041_e82_400bps_sup_v4.3.0|r1041_e82_400bps_hac_v4.2.0|r1041_e82_400bps_sup_v4.2.0|r941_sup_plant_g610|r941_min_fast_g507|r941_prom_fast_g507|r941_min_fast_g303|r941_min_high_g303|r941_min_high_g330|r941_prom_fast_g303|r941_prom_high_g303|r941_prom_high_g330|r941_min_high_g344|r941_min_high_g351|r941_min_high_g360|r941_prom_high_g344|r941_prom_high_g360|r941_prom_high_g4011|r10_min_high_g303|r10_min_high_g340|r103_min_high_g345|r103_min_high_g360|r103_prom_high_g360|r103_fast_g507|r103_hac_g507|r103_sup_g507|r104_e81_fast_g5015|r104_e81_sup_g5015|r104_e81_hac_g5015|r104_e81_sup_g610|r1041_e82_400bps_hac_g615|r1041_e82_400bps_fast_g615|r1041_e82_400bps_fast_g632|r1041_e82_260bps_fast_g632|r1041_e82_400bps_hac_g632|r1041_e82_400bps_sup_g615|r1041_e82_260bps_hac_g632|r1041_e82_260bps_sup_g632|r1041_e82_400bps_hac_v4.0.0|r1041_e82_400bps_sup_v4.0.0|r1041_e82_260bps_hac_v4.0.0|r1041_e82_260bps_sup_v4.0.0|r1041_e82_260bps_hac_v4.1.0|r1041_e82_260bps_sup_v4.1.0|r1041_e82_400bps_hac_v4.1.0|r1041_e82_400bps_sup_v4.1.0|r941_min_high_g340_rle|r941_min_hac_g507|r941_min_sup_g507|r941_prom_hac_g507|r941_prom_sup_g507|r941_e81_fast_g514|r941_e81_hac_g514|r941_e81_sup_g514]
                                   Medaka Model.  [default:
-                                  r1041_e82_400bps_sup_v4.2.0]
+                                  r1041_e82_400bps_sup_v5.0.0]
   --flyeModel [--nano-hq|--nano-corr|--nano-raw|--pacbio-raw|--pacbio-corr|--pacbio-hifi]
                                   Flye Assembly Parameter  [default: --nano-
                                   hq]
@@ -242,11 +278,16 @@ Options:
                                   be used as a database with dnaapler custom.
   --no_medaka                     Do not polish the long read assembly with
                                   Medaka.
+  --auto                          Automatically estimate the chromosome size
+                                  using KMC.
   --depth_filter FLOAT            Depth filter to pass to Plassembler. Filters
                                   out all putative plasmid contigs below this
                                   fraction of the chromosome read depth (needs
                                   to be below in both long and short read sets
                                   for hybrid).
+  --mac                           If you are running Hybracter on Mac -
+                                  installs v1.8.0 of Medaka as higher versions
+                                  break.
   --use-conda / --no-use-conda    Use conda for Snakemake rules  [default:
                                   use-conda]
   --conda-prefix PATH             Custom conda env directory
@@ -266,7 +307,7 @@ Options:
 ## `hybracter long-single`
 
 
-Run `hybracter long` on a single isolate. Instead of specifying a CSV with `--input`, use `-l` to specify the long read FASTQ file, `-s` to specify the sample name, `-c` to specify the chromosome size.
+Run `hybracter long` on a single isolate. Instead of specifying a CSV with `--input`, use `-l` to specify the long read FASTQ file, `-s` to specify the sample name, `-c` to specify the chromosome size (unless you are using `--auto`, then you don't need to specify `-c`).
 
 ```bash
 hybracter long-single -l <longread FASTQ> -s <sample name> -c <chromosome size>  -o <output_dir> -t <threads>  [other arguments]
@@ -276,24 +317,33 @@ hybracter long-single -l <longread FASTQ> -s <sample name> -c <chromosome size>
 Usage: hybracter long-single [OPTIONS] [SNAKE_ARGS]...
 
   Run hybracter long on 1 isolate
-
 Options:
   -l, --longreads TEXT            FASTQ file of longreads  [required]
   -s, --sample TEXT               Sample name.  [default: sample]
   -c, --chromosome INTEGER        FApproximate lower-bound chromosome length
                                   (in base pairs).  [default: 1000000]
   -o, --output PATH               Output directory  [default: hybracter_out]
-  --configfile TEXT               Custom config file [default:
-                                  (outputDir)/config.yaml]
+  --configfile TEXT               Custom config file [default: config.yaml]
   -t, --threads INTEGER           Number of threads to use  [default: 1]
-  --min_length INTEGER            min read length for long reads
-  --min_quality INTEGER           min read quality for long reads
-  --skip_qc                       Do not run porechop, filtlong and fastp to
-                                  QC the reads
+  --min_length INTEGER            min read length for long reads  [default:
+                                  1000]
+  --min_quality INTEGER           min read quality score for long reads in bp.
+                                  [default: 9]
+  --skip_qc                       Do not run porechop_abi, filtlong and fastp
+                                  to QC the reads
   -d, --databases PATH            Plassembler Databases directory.
-  --medakaModel [r1041_e82_400bps_hac_v4.2.0|r1041_e82_400bps_sup_v4.2.0|r941_sup_plant_g610|r941_min_fast_g507|r941_prom_fast_g507|r941_min_fast_g303|r941_min_high_g303|r941_min_high_g330|r941_prom_fast_g303|r941_prom_high_g303|r941_prom_high_g330|r941_min_high_g344|r941_min_high_g351|r941_min_high_g360|r941_prom_high_g344|r941_prom_high_g360|r941_prom_high_g4011|r10_min_high_g303|r10_min_high_g340|r103_min_high_g345|r103_min_high_g360|r103_prom_high_g360|r103_fast_g507|r103_hac_g507|r103_sup_g507|r104_e81_fast_g5015|r104_e81_sup_g5015|r104_e81_hac_g5015|r104_e81_sup_g610|r1041_e82_400bps_hac_g615|r1041_e82_400bps_fast_g615|r1041_e82_400bps_fast_g632|r1041_e82_260bps_fast_g632|r1041_e82_400bps_hac_g632|r1041_e82_400bps_sup_g615|r1041_e82_260bps_hac_g632|r1041_e82_260bps_sup_g632|r1041_e82_400bps_hac_v4.0.0|r1041_e82_400bps_sup_v4.0.0|r1041_e82_260bps_hac_v4.0.0|r1041_e82_260bps_sup_v4.0.0|r1041_e82_260bps_hac_v4.1.0|r1041_e82_260bps_sup_v4.1.0|r1041_e82_400bps_hac_v4.1.0|r1041_e82_400bps_sup_v4.1.0|r941_min_high_g340_rle|r941_min_hac_g507|r941_min_sup_g507|r941_prom_hac_g507|r941_prom_sup_g507|r941_e81_fast_g514|r941_e81_hac_g514|r941_e81_sup_g514]
+  --subsample_depth INTEGER       subsampled long read depth to subsample with
+                                  Filtlong. By default is 100x.  [default:
+                                  100]
+  --min_depth INTEGER             minimum long read depth to continue the run.
+                                  By default is 0x. Hybracter will error and
+                                  exit if a sample has less than
+                                  min_depth*chromosome_size bases of long-
+                                  reads left AFTER filtlong and porechop-ABI
+                                  steps are run.  [default: 0]
+  --medakaModel [r1041_e82_400bps_sup_v5.0.0|r1041_e82_400bps_hac_v5.0.0|r1041_e82_400bps_hac_v4.3.0|r1041_e82_400bps_sup_v4.3.0|r1041_e82_400bps_hac_v4.2.0|r1041_e82_400bps_sup_v4.2.0|r941_sup_plant_g610|r941_min_fast_g507|r941_prom_fast_g507|r941_min_fast_g303|r941_min_high_g303|r941_min_high_g330|r941_prom_fast_g303|r941_prom_high_g303|r941_prom_high_g330|r941_min_high_g344|r941_min_high_g351|r941_min_high_g360|r941_prom_high_g344|r941_prom_high_g360|r941_prom_high_g4011|r10_min_high_g303|r10_min_high_g340|r103_min_high_g345|r103_min_high_g360|r103_prom_high_g360|r103_fast_g507|r103_hac_g507|r103_sup_g507|r104_e81_fast_g5015|r104_e81_sup_g5015|r104_e81_hac_g5015|r104_e81_sup_g610|r1041_e82_400bps_hac_g615|r1041_e82_400bps_fast_g615|r1041_e82_400bps_fast_g632|r1041_e82_260bps_fast_g632|r1041_e82_400bps_hac_g632|r1041_e82_400bps_sup_g615|r1041_e82_260bps_hac_g632|r1041_e82_260bps_sup_g632|r1041_e82_400bps_hac_v4.0.0|r1041_e82_400bps_sup_v4.0.0|r1041_e82_260bps_hac_v4.0.0|r1041_e82_260bps_sup_v4.0.0|r1041_e82_260bps_hac_v4.1.0|r1041_e82_260bps_sup_v4.1.0|r1041_e82_400bps_hac_v4.1.0|r1041_e82_400bps_sup_v4.1.0|r941_min_high_g340_rle|r941_min_hac_g507|r941_min_sup_g507|r941_prom_hac_g507|r941_prom_sup_g507|r941_e81_fast_g514|r941_e81_hac_g514|r941_e81_sup_g514]
                                   Medaka Model.  [default:
-                                  r1041_e82_400bps_sup_v4.2.0]
+                                  r1041_e82_400bps_sup_v5.0.0]
   --flyeModel [--nano-hq|--nano-corr|--nano-raw|--pacbio-raw|--pacbio-corr|--pacbio-hifi]
                                   Flye Assembly Parameter  [default: --nano-
                                   hq]
@@ -305,16 +355,16 @@ Options:
                                   be used as a database with dnaapler custom.
   --no_medaka                     Do not polish the long read assembly with
                                   Medaka.
+  --auto                          Automatically estimate the chromosome size
+                                  using KMC.
   --depth_filter FLOAT            Depth filter to pass to Plassembler. Filters
                                   out all putative plasmid contigs below this
                                   fraction of the chromosome read depth (needs
                                   to be below in both long and short read sets
                                   for hybrid).
-  --logic [best|last]             Hybracter logic to select best assembly. Use
-                                  --best to pick best assembly based on ALE
-                                  (hybrid) or pyrodigal mean length (long).
-                                  Use --last to pick the last polishing round
-                                  regardless.  [default: best]
+  --mac                           If you are running Hybracter on Mac -
+                                  installs v1.8.0 of Medaka as higher versions
+                                  break.
   --use-conda / --no-use-conda    Use conda for Snakemake rules  [default:
                                   use-conda]
   --conda-prefix PATH             Custom conda env directory
@@ -322,5 +372,10 @@ Options:
                                   --rerun-incomplete, --printshellcmds,
                                   --nolock, --show-failed-logs, --conda-
                                   frontend mamba]
+  --logic [best|last]             Hybracter logic to select best assembly. Use
+                                  --best to pick best assembly based on ALE
+                                  (hybrid) or pyrodigal mean length (long).
+                                  Use --last to pick the last polishing round
+                                  regardless.  [default: best]
   -h, --help                      Show this message and exit.
 ```
diff --git a/hybracter/__main__.py b/hybracter/__main__.py
index fbf10fee..dd27f734 100644
--- a/hybracter/__main__.py
+++ b/hybracter/__main__.py
@@ -85,6 +85,14 @@ def common_options(func):
             default=100,
             show_default=True,
         ),
+        click.option(
+            "--min_depth",
+            "min_depth",
+            help="minimum long read depth to continue the run. By default is 0x.\n Hybracter will error and exit if a sample has less than min_depth*chromosome_size bases of long-reads left AFTER filtlong and porechop-ABI steps are run.",
+            type=int,
+            default=0,
+            show_default=True,
+        ),
         click.option(
             "--medakaModel",
             "medakaModel",
@@ -130,6 +138,12 @@ def common_options(func):
             is_flag=True,
             default=False,
         ),
+        click.option(
+            "--auto",
+            help="Automatically estimate the chromosome size using KMC.",
+            is_flag=True,
+            default=False,
+        ),
         click.option(
             "--depth_filter",
             help="Depth filter to pass to Plassembler. Filters out all putative plasmid contigs below this fraction of the chromosome read depth (needs to be below in both long and short read sets for hybrid).",
@@ -164,7 +178,7 @@ def common_options(func):
                 "--printshellcmds",
                 "--nolock",
                 "--show-failed-logs",
-                "--conda-frontend mamba",
+                "--conda-frontend conda",
             ],
             help="Customise Snakemake runtime args",
             show_default=True,
@@ -344,7 +358,7 @@ def print_splash():
     ),
 )
 @click.option("-i", "--input", "_input", help="Input csv", type=str, required=True)
-@click.option("--datadir", "datadir", help="Directory where FASTQs are. Will be added to the filename if in input csv if provided", type=str, default=None)
+@click.option("--datadir", "datadir", help="Directory/ies where FASTQs are. Can specify 1 directory (long and short FASTQs in the same directory) or 2 (long and short FASTQs in separate directories).\n If you specify 2, they must be separated by a comma e.g. dirlong,dirshort.\n Will be added to the filename if an input csv if provided.", type=str, default=None)
 @click.option(
     "--no_pypolca",
     help="Do not use pypolca to polish assemblies with short reads",
@@ -369,11 +383,13 @@ def hybrid(
     _input,
     datadir,
     no_medaka,
+    auto,
     no_pypolca,
     skip_qc,
     medakaModel,
     databases,
     subsample_depth,
+    min_depth,
     min_quality,
     flyeModel,
     min_length,
@@ -398,6 +414,7 @@ def hybrid(
             "min_length": min_length,
             "databases": databases,
             "subsample_depth": subsample_depth,
+            "min_depth":min_depth,
             "min_quality": min_quality,
             "no_pypolca": no_pypolca,
             "skip_qc": skip_qc,
@@ -406,6 +423,7 @@ def hybrid(
             "contaminants": contaminants,
             "dnaapler_custom_db": dnaapler_custom_db,
             "no_medaka": no_medaka,
+            "auto": auto,
             "mac": mac,
             "depth_filter": depth_filter,
             "single": False,
@@ -495,6 +513,7 @@ def hybrid_single(
     medakaModel,
     databases,
     subsample_depth,
+    min_depth,
     min_quality,
     flyeModel,
     min_length,
@@ -502,6 +521,7 @@ def hybrid_single(
     contaminants,
     dnaapler_custom_db,
     no_medaka,
+    auto,
     mac,
     logic,
     depth_filter,
@@ -523,6 +543,7 @@ def hybrid_single(
             "min_length": min_length,
             "databases": databases,
             "subsample_depth": subsample_depth,
+            "min_depth":min_depth,
             "min_quality": min_quality,
             "no_pypolca": no_pypolca,
             "skip_qc": skip_qc,
@@ -531,6 +552,7 @@ def hybrid_single(
             "contaminants": contaminants,
             "dnaapler_custom_db": dnaapler_custom_db,
             "no_medaka": no_medaka,
+            "auto": auto,
             "mac": mac,
             "depth_filter": depth_filter,
             "single": True,
@@ -583,6 +605,7 @@ def long(
     medakaModel,
     databases,
     subsample_depth,
+    min_depth,
     skip_qc,
     flyeModel,
     min_length,
@@ -591,6 +614,7 @@ def long(
     contaminants,
     dnaapler_custom_db,
     no_medaka,
+    auto,
     mac,
     logic,
     depth_filter,
@@ -611,11 +635,13 @@ def long(
             "skip_qc": skip_qc,
             "databases": databases,
             "subsample_depth": subsample_depth,
+            "min_depth": min_depth,
             "medakaModel": medakaModel,
             "flyeModel": flyeModel,
             "contaminants": contaminants,
             "dnaapler_custom_db": dnaapler_custom_db,
             "no_medaka": no_medaka,
+            "auto": auto,
             "mac": mac,
             "depth_filter": depth_filter,
             "single": False,
@@ -681,6 +707,7 @@ def long_single(
     medakaModel,
     databases,
     subsample_depth,
+    min_depth,
     skip_qc,
     flyeModel,
     min_length,
@@ -690,6 +717,7 @@ def long_single(
     dnaapler_custom_db,
     log,
     no_medaka,
+    auto,
     mac,
     logic,
     depth_filter,
@@ -710,11 +738,13 @@ def long_single(
             "skip_qc": skip_qc,
             "databases": databases,
             "subsample_depth": subsample_depth,
+            "min_depth": min_depth,
             "medakaModel": medakaModel,
             "flyeModel": flyeModel,
             "contaminants": contaminants,
             "dnaapler_custom_db": dnaapler_custom_db,
             "no_medaka": no_medaka,
+            "auto": auto,
             "mac": mac,
             "depth_filter": depth_filter,
             "single": True,
@@ -759,7 +789,7 @@ def long_single(
         "--printshellcmds",
         "--nolock",
         "--show-failed-logs",
-        "--conda-frontend mamba",
+        "--conda-frontend conda",
     ],
     help="Customise Snakemake runtime args",
     show_default=True,
@@ -857,6 +887,7 @@ def test_hybrid(
     output,
     log,
     no_medaka,
+    auto,
     min_length,
     min_quality,
     skip_qc,
@@ -864,6 +895,7 @@ def test_hybrid(
     flyeModel,
     databases,
     subsample_depth,
+    min_depth,
     no_pypolca,
     mac,
     contaminants,
@@ -886,10 +918,12 @@ def test_hybrid(
             "flyeModel": flyeModel,
             "databases": databases,
             "subsample_depth": subsample_depth,
+            "min_depth": min_depth,
             "no_pypolca": no_pypolca,
             "contaminants": contaminants,
             "dnaapler_custom_db": dnaapler_custom_db,
             "no_medaka": no_medaka,
+            "auto": auto,
             "mac": mac,
             "logic": logic,
             "depth_filter": depth_filter,
@@ -940,10 +974,12 @@ def test_long(
     medakaModel,
     databases,
     subsample_depth,
+    min_depth,
     flyeModel,
     contaminants,
     dnaapler_custom_db,
     no_medaka,
+    auto,
     mac,
     logic,
     depth_filter,
@@ -962,10 +998,12 @@ def test_long(
             "medakaModel": medakaModel,
             "databases": databases,
             "subsample_depth": subsample_depth,
+            "min_depth": min_depth,
             "flyeModel": flyeModel,
             "contaminants": contaminants,
             "dnaapler_custom_db": dnaapler_custom_db,
             "no_medaka": no_medaka,
+            "auto": auto,
             "mac": mac,
             "logic": logic,
             "depth_filter": depth_filter,
diff --git a/hybracter/hybracter.VERSION b/hybracter/hybracter.VERSION
index 8adc70fd..899f24fc 100644
--- a/hybracter/hybracter.VERSION
+++ b/hybracter/hybracter.VERSION
@@ -1,1 +1,1 @@
-0.8.0
\ No newline at end of file
+0.9.0
\ No newline at end of file
diff --git a/hybracter/workflow/envs/dnaapler.yaml b/hybracter/workflow/envs/dnaapler.yaml
index c09625d4..dc4d697b 100644
--- a/hybracter/workflow/envs/dnaapler.yaml
+++ b/hybracter/workflow/envs/dnaapler.yaml
@@ -3,5 +3,5 @@ channels:
     - conda-forge
     - bioconda
 dependencies:
-    - dnaapler >=0.7.0
+    - dnaapler >=0.8.1
 
diff --git a/hybracter/workflow/envs/kmc.yaml b/hybracter/workflow/envs/kmc.yaml
new file mode 100644
index 00000000..5021cca7
--- /dev/null
+++ b/hybracter/workflow/envs/kmc.yaml
@@ -0,0 +1,7 @@
+name: kmc
+channels:
+    - conda-forge
+    - bioconda
+dependencies: 
+    - kmc >=3.2.4
+    
diff --git a/hybracter/workflow/hybrid.smk b/hybracter/workflow/hybrid.smk
index 109f5c03..105d40d6 100644
--- a/hybracter/workflow/hybrid.smk
+++ b/hybracter/workflow/hybrid.smk
@@ -32,7 +32,6 @@ configfile: os.path.join(workflow.basedir, "../", "config", "config.yaml")
 config = ap.AttrMap(config)
 
 
-
 # directories
 CHECKDB = True  # to check db installations inside directories.smk
 
@@ -67,27 +66,33 @@ FLYE_MODEL = config.args.flyeModel
 LOGIC = config.args.logic
 DEPTH_FILTER = config.args.depth_filter
 SUBSAMPLE_DEPTH = config.args.subsample_depth
+MIN_DEPTH = config.args.min_depth
+AUTO = config.args.auto
 MAC = config.args.mac
 
 # MAC medaka
 
-new_models = ['r1041_e82_400bps_sup_v5.0.0',
-    'r1041_e82_400bps_hac_v5.0.0',
-    'r1041_e82_400bps_hac_v4.3.0',
-    'r1041_e82_400bps_sup_v4.3.0']
+new_models = [
+    "r1041_e82_400bps_sup_v5.0.0",
+    "r1041_e82_400bps_hac_v5.0.0",
+    "r1041_e82_400bps_hac_v4.3.0",
+    "r1041_e82_400bps_sup_v4.3.0",
+]
 
 if MAC:
     if MEDAKA_MODEL in new_models:
-        print(f"{MEDAKA_MODEL} is not available in medaka v1.8.0 as it is too new. If you want this model, try Hybracter on a Linux machine.")
+        print(
+            f"{MEDAKA_MODEL} is not available in medaka v1.8.0 as it is too new. If you want this model, try Hybracter on a Linux machine."
+        )
         print(f"Changing the medaka model to r1041_e82_400bps_sup_v4.2.0.")
         MEDAKA_MODEL = "r1041_e82_400bps_sup_v4.2.0"
-        
+
 
 # Parse the samples and read files
 
 # for hybracter hybrid
 if config.args.single is False:
-    dictReads = parseSamples(INPUT, False, SUBSAMPLE_DEPTH, DATADIR)  # long flag false
+    dictReads = parseSamples(INPUT, False, SUBSAMPLE_DEPTH, DATADIR, MIN_DEPTH, AUTO)  # long flag false
     SAMPLES = list(dictReads.keys())
 # for hybracter hybrid-single
 else:
@@ -96,6 +101,7 @@ else:
     dictReads[config.args.sample]["LR"] = config.args.longreads
     dictReads[config.args.sample]["MinChromLength"] = config.args.chromosome
     dictReads[config.args.sample]["TargetBases"] = SUBSAMPLE_DEPTH * config.args.chromosome
+    dictReads[config.args.sample]["MinBases"] = MIN_DEPTH * config.args.chromosome
     dictReads[config.args.sample]["R1"] = config.args.short_one
     dictReads[config.args.sample]["R2"] = config.args.short_two
     SAMPLES = [config.args.sample]
@@ -138,6 +144,8 @@ include: os.path.join("rules", "completeness", "aggregate.smk")
 # for short read polishing --careful
 include: os.path.join("rules", "processing", "coverage.smk")
 
+# kmc - needs to be included due to the checkpointing
+include: os.path.join("rules", "processing", "estimate_chromosome.smk")
 
 ### medaka vs no medaka
 # default - medaka will be run
diff --git a/hybracter/workflow/long.smk b/hybracter/workflow/long.smk
index a998eb9e..b48c8923 100644
--- a/hybracter/workflow/long.smk
+++ b/hybracter/workflow/long.smk
@@ -64,18 +64,24 @@ FLYE_MODEL = config.args.flyeModel
 LOGIC = config.args.logic
 DEPTH_FILTER = config.args.depth_filter
 SUBSAMPLE_DEPTH = config.args.subsample_depth
+MIN_DEPTH = config.args.min_depth
+AUTO = config.args.auto
 MAC = config.args.mac
 
 # MAC medaka
 
-new_models = ['r1041_e82_400bps_sup_v5.0.0',
-    'r1041_e82_400bps_hac_v5.0.0',
-    'r1041_e82_400bps_hac_v4.3.0',
-    'r1041_e82_400bps_sup_v4.3.0']
+new_models = [
+    "r1041_e82_400bps_sup_v5.0.0",
+    "r1041_e82_400bps_hac_v5.0.0",
+    "r1041_e82_400bps_hac_v4.3.0",
+    "r1041_e82_400bps_sup_v4.3.0",
+]
 
 if MAC:
     if MEDAKA_MODEL in new_models:
-        print(f"{MEDAKA_MODEL} is not available in medaka v1.8.0 as it is too new. If you want this model, try Hybracter on a Linux machine.")
+        print(
+            f"{MEDAKA_MODEL} is not available in medaka v1.8.0 as it is too new. If you want this model, try Hybracter on a Linux machine."
+        )
         print(f"Changing the medaka model to r1041_e82_400bps_sup_v4.2.0.")
         MEDAKA_MODEL = "r1041_e82_400bps_sup_v4.2.0"
 
@@ -83,7 +89,7 @@ if MAC:
 
 # for hybracter hybrid
 if config.args.single is False:
-    dictReads = parseSamples(INPUT, True, SUBSAMPLE_DEPTH, DATADIR)  # long flag true
+    dictReads = parseSamples(INPUT, True, SUBSAMPLE_DEPTH, DATADIR, MIN_DEPTH, AUTO)  # long flag true
     SAMPLES = list(dictReads.keys())
 # for hybracter hybrid-single
 else:
@@ -93,12 +99,16 @@ else:
     dictReads[config.args.sample]["MinChromLength"] = config.args.chromosome
     # add target bases for filtlong
     dictReads[config.args.sample]["TargetBases"] = SUBSAMPLE_DEPTH * config.args.chromosome
+    dictReads[config.args.sample]["MinBases"] = MIN_DEPTH * config.args.chromosome
     SAMPLES = [config.args.sample]
 
 ##############################
 # Import rules and functions
 ##############################
 
+# kmc - needs to be included due to the
+include: os.path.join("rules", "processing", "estimate_chromosome.smk")
+
 # qc and host
 # depends on whehter --contaminants has been specified and --skip_qc flag activiated
 if config.args.contaminants != "none":  # where --contaminants specified
diff --git a/hybracter/workflow/rules/assembly/assemble.smk b/hybracter/workflow/rules/assembly/assemble.smk
index fd172029..61953042 100644
--- a/hybracter/workflow/rules/assembly/assemble.smk
+++ b/hybracter/workflow/rules/assembly/assemble.smk
@@ -5,6 +5,7 @@ rule assemble:
     """
     input:
         fastq=os.path.join(dir.out.qc, "{sample}_filt_trim.fastq.gz"),
+        lr_coverage=os.path.join(dir.out.coverage, "{sample}_lr.txt"), # make sure coverage is run before assembly in case below min_depth
     output:
         fasta=os.path.join(dir.out.assemblies, "{sample}", "assembly.fasta"),
         info=os.path.join(dir.out.assemblies, "{sample}", "assembly_info.txt"),
diff --git a/hybracter/workflow/rules/assembly/plassembler.smk b/hybracter/workflow/rules/assembly/plassembler.smk
index ccae9a95..52392815 100644
--- a/hybracter/workflow/rules/assembly/plassembler.smk
+++ b/hybracter/workflow/rules/assembly/plassembler.smk
@@ -9,6 +9,7 @@ rule plassembler_hybrid:
         l=os.path.join(dir.out.qc, "{sample}_filt_trim.fastq.gz"),
         r1=os.path.join(dir.out.fastp, "{sample}_1.fastq.gz"),
         r2=os.path.join(dir.out.fastp, "{sample}_2.fastq.gz"),
+        kmc=os.path.join(dir.out.kmc,"{sample}", "{sample}_kmcLOG.txt")
     output:
         fasta=os.path.join(
             dir.out.plassembler, "{sample}", "plassembler_plasmids.fasta"
@@ -17,7 +18,7 @@ rule plassembler_hybrid:
         version=os.path.join(dir.out.versions, "{sample}", "plassembler.version"),
     params:
         db=dir.plassemblerdb,
-        chromlen=getMinChromLength,
+        chromlen = lambda wildcards: str(getMinChromLength(kmc_log_path=os.path.join(dir.out.kmc, f"{wildcards.sample}", f"{wildcards.sample}_kmcLOG.txt"), sample=wildcards.sample,auto=AUTO)),
         outdir=os.path.join(dir.out.plassembler, "{sample}"),
         flye_dir=os.path.join(dir.out.assemblies, "{sample}"),
         depth_filter=DEPTH_FILTER,
diff --git a/hybracter/workflow/rules/assembly/plassembler_long.smk b/hybracter/workflow/rules/assembly/plassembler_long.smk
index 3dfaf45d..f1819ce0 100644
--- a/hybracter/workflow/rules/assembly/plassembler_long.smk
+++ b/hybracter/workflow/rules/assembly/plassembler_long.smk
@@ -3,7 +3,8 @@ rule plassembler_long:
     runs plassembler for long
     """
     input:
-        l=os.path.join(dir.out.qc, "{sample}_filt_trim.fastq.gz"),
+        kmc=os.path.join(dir.out.kmc,"{sample}", "{sample}_kmcLOG.txt"),
+        l=os.path.join(dir.out.qc, "{sample}_filt_trim.fastq.gz")
     output:
         fasta=os.path.join(
             dir.out.plassembler, "{sample}", "plassembler_plasmids.fasta"
@@ -12,7 +13,7 @@ rule plassembler_long:
         version=os.path.join(dir.out.versions, "{sample}", "plassembler.version"),
     params:
         db=dir.plassemblerdb,
-        chromlen=getMinChromLength,
+        chromlen = lambda wildcards: str(getMinChromLength(kmc_log_path=os.path.join(dir.out.kmc, f"{wildcards.sample}", f"{wildcards.sample}_kmcLOG.txt"), sample=wildcards.sample,auto=AUTO)),
         outdir=os.path.join(dir.out.plassembler, "{sample}"),
         flye_dir=os.path.join(dir.out.assemblies, "{sample}"),
         depth_filter=DEPTH_FILTER,
diff --git a/hybracter/workflow/rules/finalise/select_best_assembly.smk b/hybracter/workflow/rules/finalise/select_best_assembly.smk
index 9c6ea1a5..95864b3e 100644
--- a/hybracter/workflow/rules/finalise/select_best_assembly.smk
+++ b/hybracter/workflow/rules/finalise/select_best_assembly.smk
@@ -102,7 +102,7 @@ rule dnaapler_pre_chrom:
         os.path.join(dir.out.stderr, "dnaapler", "{sample}_pre_chrom.log"),
     shell:
         """
-        dnaapler all -i {input.fasta} -o {params.dir} --ignore {input.ignore_list} -p {wildcards.sample} -t {threads} -a nearest --db dnaa,repa -f 2> {log}
+        dnaapler all -i {input.fasta} -o {params.dir} --ignore {input.ignore_list} -p {wildcards.sample} -t {threads} -a nearest --db dnaa,repa,cog1474 -f 2> {log}
         """
 
 
diff --git a/hybracter/workflow/rules/finalise/select_best_assembly_long.smk b/hybracter/workflow/rules/finalise/select_best_assembly_long.smk
index c7df57f6..129b8c7f 100644
--- a/hybracter/workflow/rules/finalise/select_best_assembly_long.smk
+++ b/hybracter/workflow/rules/finalise/select_best_assembly_long.smk
@@ -49,7 +49,7 @@ rule dnaapler_pre_chrom:
         os.path.join(dir.out.stderr, "dnaapler", "{sample}_pre_chrom.log"),
     shell:
         """
-        dnaapler all -i {input.fasta} -o {params.dir} --ignore {input.ignore_list} -p {wildcards.sample} -t {threads} -a nearest --db dnaa,repa -f 2> {log}
+        dnaapler all -i {input.fasta} -o {params.dir} --ignore {input.ignore_list} -p {wildcards.sample} -t {threads} -a nearest --db dnaa,repa,cog1474 -f 2> {log}
         """
 
 
diff --git a/hybracter/workflow/rules/finalise/select_best_assembly_long_no_medaka.smk b/hybracter/workflow/rules/finalise/select_best_assembly_long_no_medaka.smk
index e6777be1..dfa70ef8 100644
--- a/hybracter/workflow/rules/finalise/select_best_assembly_long_no_medaka.smk
+++ b/hybracter/workflow/rules/finalise/select_best_assembly_long_no_medaka.smk
@@ -30,7 +30,7 @@ rule dnaapler_pre_chrom:
         os.path.join(dir.out.stderr, "dnaapler", "{sample}_pre_chrom.log"),
     shell:
         """
-        dnaapler all -i {input.fasta} -o {params.dir} --ignore {input.ignore_list} -p {wildcards.sample} -t {threads} -a nearest --db dnaa,repa -f 2> {log}
+        dnaapler all -i {input.fasta} -o {params.dir} --ignore {input.ignore_list} -p {wildcards.sample} -t {threads} -a nearest --db dnaa,repa,cog1474 -f 2> {log}
         """
 
 
diff --git a/hybracter/workflow/rules/finalise/select_best_assembly_no_medaka.smk b/hybracter/workflow/rules/finalise/select_best_assembly_no_medaka.smk
index db84708d..731fc528 100644
--- a/hybracter/workflow/rules/finalise/select_best_assembly_no_medaka.smk
+++ b/hybracter/workflow/rules/finalise/select_best_assembly_no_medaka.smk
@@ -98,7 +98,7 @@ rule dnaapler_pre_chrom:
         os.path.join(dir.out.stderr, "dnaapler", "{sample}_pre_chrom.log"),
     shell:
         """
-        dnaapler all -i {input.fasta} -o {params.dir} --ignore {input.ignore_list} -p {wildcards.sample} -t {threads} -a nearest --db dnaa,repa -f 2> {log}
+        dnaapler all -i {input.fasta} -o {params.dir} --ignore {input.ignore_list} -p {wildcards.sample} -t {threads} -a nearest --db dnaa,repa,cog1474 -f 2> {log}
         """
 
 
diff --git a/hybracter/workflow/rules/polishing/long_read_polish.smk b/hybracter/workflow/rules/polishing/long_read_polish.smk
index 8d9298c6..52231952 100644
--- a/hybracter/workflow/rules/polishing/long_read_polish.smk
+++ b/hybracter/workflow/rules/polishing/long_read_polish.smk
@@ -50,13 +50,14 @@ rule medaka_round_1_extract_intermediate_assembly:
         fasta=os.path.join(dir.out.medaka_rd_1, "{sample}", "consensus.fasta"),
         completeness_check=os.path.join(dir.out.completeness, "{sample}.txt"),
         info=os.path.join(dir.out.assemblies, "{sample}", "assembly_info.txt"),
+        kmc=os.path.join(dir.out.kmc,"{sample}", "{sample}_kmcLOG.txt")
     output:
         fasta=os.path.join(
             dir.out.intermediate_assemblies, "{sample}", "{sample}_medaka_rd_1.fasta"
         ),
         ignore_list=os.path.join(dir.out.medaka_rd_1, "{sample}_ignore_list.txt"),
     params:
-        min_chrom_length=getMinChromLength,
+        min_chrom_length=lambda wildcards: str(getMinChromLength(kmc_log_path=os.path.join(dir.out.kmc, f"{wildcards.sample}", f"{wildcards.sample}_kmcLOG.txt"), sample=wildcards.sample,auto=AUTO)),
         polypolish_flag=False,
     conda:
         os.path.join(dir.env, "scripts.yaml")
@@ -186,7 +187,7 @@ rule medaka_round_2_extract_intermediate_assembly:
         ),
         ignore_list=os.path.join(dir.out.medaka_rd_2, "{sample}_ignore_list.txt"),
     params:
-        min_chrom_length=getMinChromLength,
+        min_chrom_length=lambda wildcards: str(getMinChromLength(kmc_log_path=os.path.join(dir.out.kmc, f"{wildcards.sample}", f"{wildcards.sample}_kmcLOG.txt"), sample=wildcards.sample,auto=AUTO)),
         polypolish_flag=False,
     conda:
         os.path.join(dir.env, "scripts.yaml")
diff --git a/hybracter/workflow/rules/polishing/short_read_polish.smk b/hybracter/workflow/rules/polishing/short_read_polish.smk
index 742e9a3d..726abe01 100644
--- a/hybracter/workflow/rules/polishing/short_read_polish.smk
+++ b/hybracter/workflow/rules/polishing/short_read_polish.smk
@@ -53,7 +53,7 @@ rule polypolish:
         fasta=os.path.join(dir.out.medaka_rd_2, "{sample}", "consensus.fasta"),
         sam1=os.path.join(dir.out.bwa, "{sample}_1.sam"),
         sam2=os.path.join(dir.out.bwa, "{sample}_2.sam"),
-        coverage=os.path.join(dir.out.coverage, "{sample}.txt"),
+        coverage=os.path.join(dir.out.coverage, "{sample}_sr.txt"),
     output:
         fasta=os.path.join(dir.out.polypolish, "{sample}.fasta"),
         version=os.path.join(dir.out.versions, "{sample}", "polypolish.version"),
@@ -88,13 +88,14 @@ rule polypolish_extract_intermediate_assembly:
         fasta=os.path.join(dir.out.polypolish, "{sample}.fasta"),
         completeness_check=os.path.join(dir.out.completeness, "{sample}.txt"),
         info=os.path.join(dir.out.assemblies, "{sample}", "assembly_info.txt"),
+        kmc=os.path.join(dir.out.kmc,"{sample}", "{sample}_kmcLOG.txt")
     output:
         fasta=os.path.join(
             dir.out.intermediate_assemblies, "{sample}", "{sample}_polypolish.fasta"
         ),
         ignore_list=os.path.join(dir.out.polypolish, "{sample}_ignore_list.txt"),
     params:
-        min_chrom_length=getMinChromLength,
+        min_chrom_length=lambda wildcards: str(getMinChromLength(kmc_log_path=os.path.join(dir.out.kmc, f"{wildcards.sample}", f"{wildcards.sample}_kmcLOG.txt"), sample=wildcards.sample,auto=AUTO)),
         polypolish_flag=True,
     conda:
         os.path.join(dir.env, "scripts.yaml")
diff --git a/hybracter/workflow/rules/polishing/short_read_polish_incomplete.smk b/hybracter/workflow/rules/polishing/short_read_polish_incomplete.smk
index a8f3a7e2..82503315 100644
--- a/hybracter/workflow/rules/polishing/short_read_polish_incomplete.smk
+++ b/hybracter/workflow/rules/polishing/short_read_polish_incomplete.smk
@@ -53,7 +53,7 @@ rule polypolish_incomplete:
         fasta=os.path.join(dir.out.medaka_incomplete, "{sample}", "consensus.fasta"),
         sam1=os.path.join(dir.out.bwa_incomplete, "{sample}_1.sam"),
         sam2=os.path.join(dir.out.bwa_incomplete, "{sample}_2.sam"),
-        coverage=os.path.join(dir.out.coverage, "{sample}.txt"),
+        coverage=os.path.join(dir.out.coverage, "{sample}_sr.txt"),
     output:
         fasta=os.path.join(dir.out.polypolish_incomplete, "{sample}.fasta"),
         version=os.path.join(
diff --git a/hybracter/workflow/rules/polishing/short_read_polish_incomplete_no_medaka.smk b/hybracter/workflow/rules/polishing/short_read_polish_incomplete_no_medaka.smk
index e857c92e..ee9cc51b 100644
--- a/hybracter/workflow/rules/polishing/short_read_polish_incomplete_no_medaka.smk
+++ b/hybracter/workflow/rules/polishing/short_read_polish_incomplete_no_medaka.smk
@@ -53,7 +53,7 @@ rule polypolish_incomplete:
         fasta=os.path.join(dir.out.incomp_pre_polish, "{sample}.fasta"),
         sam1=os.path.join(dir.out.bwa_incomplete, "{sample}_1.sam"),
         sam2=os.path.join(dir.out.bwa_incomplete, "{sample}_2.sam"),
-        coverage=os.path.join(dir.out.coverage, "{sample}.txt"),
+        coverage=os.path.join(dir.out.coverage, "{sample}_sr.txt"),
     output:
         fasta=os.path.join(dir.out.polypolish_incomplete, "{sample}.fasta"),
         version=os.path.join(
diff --git a/hybracter/workflow/rules/polishing/short_read_polish_no_medaka.smk b/hybracter/workflow/rules/polishing/short_read_polish_no_medaka.smk
index 0f12df0e..74486570 100644
--- a/hybracter/workflow/rules/polishing/short_read_polish_no_medaka.smk
+++ b/hybracter/workflow/rules/polishing/short_read_polish_no_medaka.smk
@@ -73,7 +73,7 @@ rule polypolish:
         ),
         sam1=os.path.join(dir.out.bwa, "{sample}_1.sam"),
         sam2=os.path.join(dir.out.bwa, "{sample}_2.sam"),
-        coverage=os.path.join(dir.out.coverage, "{sample}.txt"),
+        coverage=os.path.join(dir.out.coverage, "{sample}_sr.txt"),
     output:
         fasta=os.path.join(dir.out.polypolish, "{sample}.fasta"),
         version=os.path.join(dir.out.versions, "{sample}", "polypolish.version"),
@@ -108,13 +108,14 @@ rule polypolish_extract_intermediate_assembly:
         fasta=os.path.join(dir.out.polypolish, "{sample}.fasta"),
         completeness_check=os.path.join(dir.out.completeness, "{sample}.txt"),
         info=os.path.join(dir.out.assemblies, "{sample}", "assembly_info.txt"),
+        kmc=os.path.join(dir.out.kmc,"{sample}", "{sample}_kmcLOG.txt")
     output:
         fasta=os.path.join(
             dir.out.intermediate_assemblies, "{sample}", "{sample}_polypolish.fasta"
         ),
         ignore_list=os.path.join(dir.out.polypolish, "{sample}_ignore_list.txt"),
     params:
-        min_chrom_length=getMinChromLength,
+        min_chrom_length=lambda wildcards: str(getMinChromLength(kmc_log_path=os.path.join(dir.out.kmc, f"{wildcards.sample}", f"{wildcards.sample}_kmcLOG.txt"), sample=wildcards.sample,auto=AUTO)),
         polypolish_flag=True,
     conda:
         os.path.join(dir.env, "scripts.yaml")
diff --git a/hybracter/workflow/rules/polishing/short_read_pypolca.smk b/hybracter/workflow/rules/polishing/short_read_pypolca.smk
index 66cd5432..8b700103 100644
--- a/hybracter/workflow/rules/polishing/short_read_pypolca.smk
+++ b/hybracter/workflow/rules/polishing/short_read_pypolca.smk
@@ -8,7 +8,7 @@ rule pypolca:
         polypolish_fasta=os.path.join(dir.out.polypolish, "{sample}.fasta"),
         r1=os.path.join(dir.out.fastp, "{sample}_1.fastq.gz"),
         r2=os.path.join(dir.out.fastp, "{sample}_2.fastq.gz"),
-        coverage=os.path.join(dir.out.coverage, "{sample}.txt"),
+        coverage=os.path.join(dir.out.coverage, "{sample}_sr.txt"),
     output:
         fasta=os.path.join(dir.out.pypolca, "{sample}", "{sample}_corrected.fasta"),
         version=os.path.join(dir.out.versions, "{sample}", "pypolca_complete.version"),
@@ -47,13 +47,14 @@ rule pypolca_extract_intermediate_assembly:
         fasta=os.path.join(dir.out.pypolca, "{sample}", "{sample}_corrected.fasta"),
         completeness_check=os.path.join(dir.out.completeness, "{sample}.txt"),
         info=os.path.join(dir.out.assemblies, "{sample}", "assembly_info.txt"),
+        kmc=os.path.join(dir.out.kmc,"{sample}", "{sample}_kmcLOG.txt")
     output:
         fasta=os.path.join(
             dir.out.intermediate_assemblies, "{sample}", "{sample}_pypolca.fasta"
         ),
         ignore_list=os.path.join(dir.out.pypolca, "{sample}_ignore_list.txt"),
     params:
-        min_chrom_length=getMinChromLength,
+        min_chrom_length=lambda wildcards: str(getMinChromLength(kmc_log_path=os.path.join(dir.out.kmc, f"{wildcards.sample}", f"{wildcards.sample}_kmcLOG.txt"), sample=wildcards.sample,auto=AUTO)),
         polypolish_flag=True,
     conda:
         os.path.join(dir.env, "scripts.yaml")
@@ -105,7 +106,7 @@ rule pypolca_incomplete:
         polypolish_fasta=os.path.join(dir.out.polypolish_incomplete, "{sample}.fasta"),
         r1=os.path.join(dir.out.fastp, "{sample}_1.fastq.gz"),
         r2=os.path.join(dir.out.fastp, "{sample}_2.fastq.gz"),
-        coverage=os.path.join(dir.out.coverage, "{sample}.txt"),
+        coverage=os.path.join(dir.out.coverage, "{sample}_sr.txt"),
     output:
         fasta=os.path.join(
             dir.out.pypolca_incomplete, "{sample}", "{sample}_corrected.fasta"
diff --git a/hybracter/workflow/rules/polishing/short_read_pypolca_no_medaka.smk b/hybracter/workflow/rules/polishing/short_read_pypolca_no_medaka.smk
index b3836cc8..44dd8c5f 100644
--- a/hybracter/workflow/rules/polishing/short_read_pypolca_no_medaka.smk
+++ b/hybracter/workflow/rules/polishing/short_read_pypolca_no_medaka.smk
@@ -8,7 +8,7 @@ rule pypolca:
         polypolish_fasta=os.path.join(dir.out.polypolish, "{sample}.fasta"),
         r1=os.path.join(dir.out.fastp, "{sample}_1.fastq.gz"),
         r2=os.path.join(dir.out.fastp, "{sample}_2.fastq.gz"),
-        coverage=os.path.join(dir.out.coverage, "{sample}.txt"),
+        coverage=os.path.join(dir.out.coverage, "{sample}_sr.txt"),
     output:
         fasta=os.path.join(dir.out.pypolca, "{sample}", "{sample}_corrected.fasta"),
         version=os.path.join(dir.out.versions, "{sample}", "pypolca_complete.version"),
@@ -47,13 +47,14 @@ rule pypolca_extract_intermediate_assembly:
         fasta=os.path.join(dir.out.pypolca, "{sample}", "{sample}_corrected.fasta"),
         completeness_check=os.path.join(dir.out.completeness, "{sample}.txt"),
         info=os.path.join(dir.out.assemblies, "{sample}", "assembly_info.txt"),
+        kmc=os.path.join(dir.out.kmc,"{sample}", "{sample}_kmcLOG.txt")
     output:
         fasta=os.path.join(
             dir.out.intermediate_assemblies, "{sample}", "{sample}_pypolca.fasta"
         ),
         ignore_list=os.path.join(dir.out.pypolca, "{sample}_ignore_list.txt"),
     params:
-        min_chrom_length=getMinChromLength,
+        min_chrom_length=lambda wildcards: str(getMinChromLength(kmc_log_path=os.path.join(dir.out.kmc, f"{wildcards.sample}", f"{wildcards.sample}_kmcLOG.txt"), sample=wildcards.sample,auto=AUTO)),
         polypolish_flag=True,
     conda:
         os.path.join(dir.env, "scripts.yaml")
@@ -102,7 +103,7 @@ rule pypolca_incomplete:
         polypolish_fasta=os.path.join(dir.out.polypolish_incomplete, "{sample}.fasta"),
         r1=os.path.join(dir.out.fastp, "{sample}_1.fastq.gz"),
         r2=os.path.join(dir.out.fastp, "{sample}_2.fastq.gz"),
-        coverage=os.path.join(dir.out.coverage, "{sample}.txt"),
+        coverage=os.path.join(dir.out.coverage, "{sample}_sr.txt"),
     output:
         fasta=os.path.join(
             dir.out.pypolca_incomplete, "{sample}", "{sample}_corrected.fasta"
diff --git a/hybracter/workflow/rules/preflight/directories.smk b/hybracter/workflow/rules/preflight/directories.smk
index 2ba31fb5..d7ef0241 100644
--- a/hybracter/workflow/rules/preflight/directories.smk
+++ b/hybracter/workflow/rules/preflight/directories.smk
@@ -40,6 +40,9 @@ dir.out.processing = os.path.join(dir.out.base, "processing")
 dir.out.complete = os.path.join(dir.out.processing, "complete")
 dir.out.incomplete = os.path.join(dir.out.processing, "incomplete")
 dir.out.qc = os.path.join(dir.out.processing, "qc")
+dir.out.kmc = os.path.join(dir.out.processing, "kmc")
+# auto only - write chrom size to file
+dir.out.chrom_size = os.path.join(dir.out.processing, "chrom_size")
 
 # logs and benchmarks
 dir.out.bench = os.path.join(dir.out.base, "benchmarks")
diff --git a/hybracter/workflow/rules/preflight/functions.smk b/hybracter/workflow/rules/preflight/functions.smk
index 0d4f9997..b9437039 100644
--- a/hybracter/workflow/rules/preflight/functions.smk
+++ b/hybracter/workflow/rules/preflight/functions.smk
@@ -6,6 +6,42 @@ import gzip
 from Bio import SeqIO
 import sys
 import os
+import re
+
+def getMinChromLength(kmc_log_path, sample, auto):
+    if auto:
+        with open(kmc_log_path, "r") as file:
+            for line in file:
+                if "No. of unique counted k-mers" in line:
+                    # keep 80% of kmers as lower bound for chromosome
+                    return int(float(re.search(r"No. of unique counted k-mers\s*:\s*([\d\.eE+-]+)", line).group(1)) * 0.8)
+    else:
+        return dictReads[sample]["MinChromLength"]
+
+# get min_depth 
+def getMinBases(kmc_log_path, sample, auto, min_depth):
+    if auto:
+        with open(kmc_log_path, "r") as file:
+            for line in file:
+                if "No. of unique counted k-mers" in line:
+                    # keep 80% of kmers as lower bound for chromosome
+                    chrom_size = int(float(re.search(r"No. of unique counted k-mers\s*:\s*([\d\.eE+-]+)", line).group(1)) * 0.8)
+                    return int(min_depth*chrom_size)
+    else:
+        return dictReads[sample]["MinBases"]
+
+def getTargetBases(kmc_log_path, sample, auto, target_depth):
+    if auto:
+        with open(kmc_log_path, "r") as file:
+            for line in file:
+                if "No. of unique counted k-mers" in line:
+                    # keep 80% of kmers as lower bound for chromosome
+                    chrom_size = int(float(re.search(r"No. of unique counted k-mers\s*:\s*([\d\.eE+-]+)", line).group(1)) * 0.8)
+                    return int(target_depth*chrom_size)
+    else:
+        return dictReads[sample]["TargetBases"]
+
+
 
 
 # define functions
@@ -14,16 +50,6 @@ def get_input_lr_fastqs(wildcards):
     return dictReads[wildcards.sample]["LR"]
 
 
-# get min chrom length (define chrom size)
-def getMinChromLength(wildcards):
-    return dictReads[wildcards.sample]["MinChromLength"]
-
-
-# get target bases for filtlong subsampling
-def getTargetBases(wildcards):
-    return dictReads[wildcards.sample]["TargetBases"]
-
-
 def get_input_r1(wildcards):
     return dictReads[wildcards.sample]["R1"]
 
diff --git a/hybracter/workflow/rules/preflight/samples.smk b/hybracter/workflow/rules/preflight/samples.smk
index 2adb8aa2..9ffaf0c0 100644
--- a/hybracter/workflow/rules/preflight/samples.smk
+++ b/hybracter/workflow/rules/preflight/samples.smk
@@ -5,52 +5,111 @@ Function for parsing the config and identifying samples and read files
 from itertools import chain
 import sys
 
+"""
+check_datadir function
+"""
+
+
+def check_datadir(datadir: str):
+    # if only one (hybracter long)
+    if "," not in datadir:
+        return datadir, None
+    else:
+        # if more than one (hybracter hybrid)
+        parts = datadir.split(",")
+
+        if len(parts) > 2:
+            sys.exit(
+                f"Error: You must specify a maximum of two comma separated input directories with --datadir, the first being where the long-read FASTQs are, the second the short-read.\n You have specified --datadir {datadir}.\n Please check this."
+            )
+
+        # If exactly one comma, assign the two parts
+        datadirlong, datadirshort = parts
+        return datadirlong, datadirshort
+
+
 """
 long
 """
 
 
-def samplesFromCsvLong(csvFile, subsample_depth, datadir):
+def samplesFromCsvLong(csvFile, subsample_depth, datadir, min_depth, auto):
     """
     Read samples and files from a CSV Long Read Only
     3 cols
     1 = sample
     2 = Long read
     3 = MinChromLength
+    IF auto - only need 2
+
+    chrom_size, target_bases and min_depth handled with functions if --auto
     """
     outDict = {}
     with open(csvFile, "r") as csv:
         for line in csv:
             l = line.strip().split(",")
-            if len(l) == 3:
-                outDict[l[0]] = {}
-                # where datafir isn't specified
-                if datadir is None:
-                    long_fastq = l[1]
-                else:
-                    long_fastq = os.path.join(datadir, l[1])
-                if os.path.isfile(long_fastq) and l[2].isnumeric():
-                    outDict[l[0]]["LR"] = long_fastq
-                    outDict[l[0]]["MinChromLength"] = l[2]
-                    outDict[l[0]]["TargetBases"] = int(l[2]) * subsample_depth
+            if auto:
+                if len(l) == 2:
+                    outDict[l[0]] = {}
+                    # where datafir isn't specified
+                    if datadir is None:
+                        long_fastq = l[1]
+                    else:
+                        datadirlong, datadirshort = check_datadir(datadir)
+                        long_fastq = os.path.join(datadirlong, l[1])
+                    if os.path.isfile(long_fastq):
+                        outDict[l[0]]["LR"] = long_fastq
+                    else:
+                        sys.stderr.write(
+                            "\n"
+                            f"    FATAL: Error parsing {csvFile}. {long_fastq} \n"
+                            f"    does not exist. \n"
+                            "    Check formatting, and that \n"
+                            "    file names and file paths are correct.\n"
+                            "\n"
+                        )
+                        sys.exit(1)
                 else:
                     sys.stderr.write(
                         "\n"
-                        f"    FATAL: Error parsing {csvFile}. {long_fastq} \n"
-                        f"    does not exist or  {l[2]} is not an integer. \n"
-                        "    Check formatting, and that \n"
-                        "    file names and file paths are correct.\n"
-                        "\n"
+                        f"    FATAL: Error parsing {csvFile}. Line {l} \n"
+                        f"    does not have 2 columns. \n"
+                        f"    Please check the formatting of {csvFile}. \n"
                     )
                     sys.exit(1)
+
             else:
-                sys.stderr.write(
-                    "\n"
-                    f"    FATAL: Error parsing {csvFile}. Line {l} \n"
-                    f"    does not have 3 columns. \n"
-                    f"    Please check the formatting of {csvFile}. \n"
-                )
-                sys.exit(1)
+                if len(l) == 3:
+                    outDict[l[0]] = {}
+                    # where datafir isn't specified
+                    if datadir is None:
+                        long_fastq = l[1]
+                    else:
+                        datadirlong, datadirshort = check_datadir(datadir)
+                        long_fastq = os.path.join(datadirlong, l[1])
+                    if os.path.isfile(long_fastq) and l[2].isnumeric():
+                        outDict[l[0]]["LR"] = long_fastq
+                        outDict[l[0]]["MinChromLength"] = l[2]
+                        outDict[l[0]]["TargetBases"] = int(l[2]) * subsample_depth
+                        outDict[l[0]]["MinBases"] = int(l[2]) * min_depth
+                    else:
+                        sys.stderr.write(
+                            "\n"
+                            f"    FATAL: Error parsing {csvFile}. {long_fastq} \n"
+                            f"    does not exist or  {l[2]} is not an integer. \n"
+                            "    Check formatting, and that \n"
+                            "    file names and file paths are correct.\n"
+                            "\n"
+                        )
+                        sys.exit(1)
+                else:
+                    sys.stderr.write(
+                        "\n"
+                        f"    FATAL: Error parsing {csvFile}. Line {l} \n"
+                        f"    does not have 3 columns. \n"
+                        f"    Please check the formatting of {csvFile}. \n"
+                    )
+                    sys.exit(1)
     return outDict
 
 
@@ -59,7 +118,7 @@ short
 """
 
 
-def samplesFromCsvShort(csvFile, subsample_depth, datadir):
+def samplesFromCsvShort(csvFile, subsample_depth, datadir, min_depth, auto):
     """
     Read samples and files from a CSV Hybrid
     5 cols
@@ -73,58 +132,119 @@ def samplesFromCsvShort(csvFile, subsample_depth, datadir):
     with open(csvFile, "r") as csv:
         for line in csv:
             l = line.strip().split(",")
-            if len(l) == 5:
-                outDict[l[0]] = {}
-                # where datafir isn't specified
-                if datadir is None:
-                    long_fastq = l[1]
-                    r1_fastq = l[3]
-                    r2_fastq = l[4]
-                else:
-                    long_fastq = os.path.join(datadir, l[1])
-                    r1_fastq = os.path.join(datadir, l[3])
-                    r2_fastq = os.path.join(datadir, l[4])
-
-                if (
-                    os.path.isfile(long_fastq)
-                    and l[2].isnumeric()
-                    and os.path.isfile(r1_fastq)
-                    and os.path.isfile(r2_fastq)
-                ):
-                    outDict[l[0]]["LR"] = long_fastq
-                    outDict[l[0]]["MinChromLength"] = l[2]
-                    outDict[l[0]]["R1"] = r1_fastq
-                    outDict[l[0]]["R2"] = r2_fastq
-                    outDict[l[0]]["TargetBases"] = int(l[2]) * subsample_depth
+            if auto:
+                if len(l) == 4:
+                    outDict[l[0]] = {}
+                    # where datafir isn't specified
+                    if datadir is None:
+                        long_fastq = l[1]
+                        r1_fastq = l[2]
+                        r2_fastq = l[3]
+                    else:
+                        # if the user supplies a datadir
+                        datadirlong, datadirshort = check_datadir(datadir)
+                        long_fastq = os.path.join(datadirlong, l[1])
+                        # all fastqs in 1 dir
+                        if datadirshort is None:
+                            r1_fastq = os.path.join(datadirlong, l[2])
+                            r2_fastq = os.path.join(datadirlong, l[3])
+                        # separate dirs
+                        else:
+                            r1_fastq = os.path.join(datadirshort, l[2])
+                            r2_fastq = os.path.join(datadirshort, l[3])                        
+
+                    if (
+                        os.path.isfile(long_fastq)
+                        and os.path.isfile(r1_fastq)
+                        and os.path.isfile(r2_fastq)
+                    ):
+                        outDict[l[0]]["LR"] = long_fastq
+                        outDict[l[0]]["R1"] = r1_fastq
+                        outDict[l[0]]["R2"] = r2_fastq
+                    else:
+                        sys.stderr.write(
+                            "\n"
+                            f"    FATAL: Error parsing {csvFile}. One of \n"
+                            f"    {long_fastq} or \n"
+                            f"    {r1_fastq} or \n"
+                            f"    {r2_fastq} \n"
+                            f"    does not exist. \n"
+                            "    Check formatting, and that \n"
+                            "    file names and file paths are correct.\n"
+                            "\n"
+                        )
+                        sys.exit(1)
                 else:
                     sys.stderr.write(
                         "\n"
-                        f"    FATAL: Error parsing {csvFile}. One of \n"
-                        f"    {long_fastq} or \n"
-                        f"    {r1_fastq} or \n"
-                        f"    {r2_fastq} \n"
-                        f"    does not exist or  {l[2]} is not an integer. \n"
-                        "    Check formatting, and that \n"
-                        "    file names and file paths are correct.\n"
-                        "\n"
+                        f"    FATAL: Error parsing {csvFile}. Line {l} \n"
+                        f"    does not have 5 columns. \n"
+                        f"    Please check the formatting of {csvFile}. \n"
                     )
                     sys.exit(1)
+
             else:
-                sys.stderr.write(
-                    "\n"
-                    f"    FATAL: Error parsing {csvFile}. Line {l} \n"
-                    f"    does not have 5 columns. \n"
-                    f"    Please check the formatting of {csvFile}. \n"
-                )
-                sys.exit(1)
+                if len(l) == 5:
+                    outDict[l[0]] = {}
+                    # where datafir isn't specified
+                    if datadir is None:
+                        long_fastq = l[1]
+                        r1_fastq = l[3]
+                        r2_fastq = l[4]
+                    else:
+                        # if the user supplies a datadir
+                        datadirlong, datadirshort = check_datadir(datadir)
+                        long_fastq = os.path.join(datadirlong, l[1])
+                        # all fastqs in 1 dir
+                        if datadirshort is None:
+                            r1_fastq = os.path.join(datadirlong, l[3])
+                            r2_fastq = os.path.join(datadirlong, l[4])
+                        # separate dirs
+                        else:
+                            r1_fastq = os.path.join(datadirshort, l[3])
+                            r2_fastq = os.path.join(datadirshort, l[4])                        
+                    if (
+                        os.path.isfile(long_fastq)
+                        and l[2].isnumeric()
+                        and os.path.isfile(r1_fastq)
+                        and os.path.isfile(r2_fastq)
+                    ):
+                        outDict[l[0]]["LR"] = long_fastq
+                        outDict[l[0]]["MinChromLength"] = l[2]
+                        outDict[l[0]]["R1"] = r1_fastq
+                        outDict[l[0]]["R2"] = r2_fastq
+                        outDict[l[0]]["TargetBases"] = int(l[2]) * subsample_depth
+                        outDict[l[0]]["MinBases"] = int(l[2]) * min_depth
+                    else:
+                        sys.stderr.write(
+                            "\n"
+                            f"    FATAL: Error parsing {csvFile}. One of \n"
+                            f"    {long_fastq} or \n"
+                            f"    {r1_fastq} or \n"
+                            f"    {r2_fastq} \n"
+                            f"    does not exist or  {l[2]} is not an integer. \n"
+                            "    Check formatting, and that \n"
+                            "    file names and file paths are correct.\n"
+                            "\n"
+                        )
+                        sys.exit(1)
+                else:
+                    sys.stderr.write(
+                        "\n"
+                        f"    FATAL: Error parsing {csvFile}. Line {l} \n"
+                        f"    does not have 5 columns. \n"
+                        f"    Please check the formatting of {csvFile}. \n"
+                    )
+                    sys.exit(1)
+
     return outDict
 
 
-def parseSamples(csvfile, long_flag, subsample_depth, datadir):
+def parseSamples(csvfile, long_flag, subsample_depth, datadir, min_depth, auto):
     if os.path.isfile(csvfile) and long_flag is True:
-        sampleDict = samplesFromCsvLong(csvfile, subsample_depth, datadir)
+        sampleDict = samplesFromCsvLong(csvfile, subsample_depth, datadir, min_depth, auto)
     elif os.path.isfile(csvfile) and long_flag is False:
-        sampleDict = samplesFromCsvShort(csvfile, subsample_depth, datadir)
+        sampleDict = samplesFromCsvShort(csvfile, subsample_depth, datadir, min_depth, auto)
     else:
         sys.stderr.write(
             "\n"
diff --git a/hybracter/workflow/rules/preflight/targets_hybrid.smk b/hybracter/workflow/rules/preflight/targets_hybrid.smk
index 5a4716d1..14dcbd99 100644
--- a/hybracter/workflow/rules/preflight/targets_hybrid.smk
+++ b/hybracter/workflow/rules/preflight/targets_hybrid.smk
@@ -23,6 +23,7 @@ hybrid
 TargetFilesHybrid = [
     os.path.join(dir.out.flags, "aggr_long_qc.flag"),
     os.path.join(dir.out.flags, "aggr_short_qc.flag"),
+    os.path.join(dir.out.flags, "aggr_kmc.flag"),
     os.path.join(dir.out.flags, "aggr_seqkit_short.flag"),
     os.path.join(dir.out.flags, "aggr_seqkit_long.flag"),
     os.path.join(dir.out.flags, "aggr_assemble.flag"),
diff --git a/hybracter/workflow/rules/preflight/targets_long.smk b/hybracter/workflow/rules/preflight/targets_long.smk
index 20a399a2..9f6835c1 100644
--- a/hybracter/workflow/rules/preflight/targets_long.smk
+++ b/hybracter/workflow/rules/preflight/targets_long.smk
@@ -9,6 +9,7 @@ long
 TargetFilesLong = [
     os.path.join(dir.out.flags, "aggr_long_qc.flag"),
     os.path.join(dir.out.flags, "aggr_seqkit_long.flag"),
+    os.path.join(dir.out.flags, "aggr_kmc.flag"),
     os.path.join(dir.out.flags, "aggr_assemble.flag"),
     os.path.join(dir.out.flags, "aggr_long_read_polish.flag"),
     os.path.join(dir.out.flags, "aggr_plassembler.flag"),
diff --git a/hybracter/workflow/rules/processing/coverage.smk b/hybracter/workflow/rules/processing/coverage.smk
index 29be969e..b24042ab 100644
--- a/hybracter/workflow/rules/processing/coverage.smk
+++ b/hybracter/workflow/rules/processing/coverage.smk
@@ -65,12 +65,37 @@ rule aggr_seqkit_short:
         """
 
 
+rule estimate_lr_coverage:
+    """
+    quickly estimates long read coverage for min_depth
+    """
+    input:
+        long_bases=os.path.join(dir.out.seqkit, "{sample}_long.txt"),
+        kmc=os.path.join(dir.out.kmc,"{sample}", "{sample}_kmcLOG.txt")
+    output:
+        lr_coverage=os.path.join(dir.out.coverage, "{sample}_lr.txt"),
+    params:
+        min_chrom_length=lambda wildcards: str(getMinChromLength(kmc_log_path=os.path.join(dir.out.kmc, f"{wildcards.sample}", f"{wildcards.sample}_kmcLOG.txt"), sample=wildcards.sample,auto=AUTO)),
+        min_bases=lambda wildcards: str(getMinBases(kmc_log_path=os.path.join(dir.out.kmc,f"{wildcards.sample}", f"{wildcards.sample}_kmcLOG.txt"),sample=wildcards.sample, auto=AUTO, min_depth=MIN_DEPTH)),
+    conda:
+        os.path.join(dir.env, "scripts.yaml")
+    resources:
+        mem_mb=config.resources.sml.mem,
+        mem=str(config.resources.sml.mem) + "MB",
+        time=config.resources.sml.time,
+    threads: config.resources.sml.cpu
+    script:
+        os.path.join(dir.scripts, "lr_coverage.py")
+        
+
+
 rule aggr_seqkit_long:
     """
-    aggregates the seqkit stats over all samples
+    aggregates the seqkit stats and coverage check over all samples
     """
     input:
         expand(os.path.join(dir.out.seqkit, "{sample}_long.txt"), sample=SAMPLES),
+        expand(os.path.join(dir.out.coverage, "{sample}_lr.txt"), sample=SAMPLES)
     output:
         flag=os.path.join(dir.out.flags, "aggr_seqkit_long.flag"),
     resources:
@@ -83,18 +108,18 @@ rule aggr_seqkit_long:
         touch {output.flag}
         """
 
-
 rule estimate_sr_coverage:
     """
-    quickly estimates short read coverage for careful
+    quickly estimates short read coverage for pypolca and polypolish careful
     """
     input:
         r1=os.path.join(dir.out.seqkit, "{sample}_r1.txt"),
         r2=os.path.join(dir.out.seqkit, "{sample}_r2.txt"),
+        kmc=os.path.join(dir.out.kmc,"{sample}", "{sample}_kmcLOG.txt")
     output:
-        sr_coverage=os.path.join(dir.out.coverage, "{sample}.txt"),
+        sr_coverage=os.path.join(dir.out.coverage, "{sample}_sr.txt"),
     params:
-        chromlen=getMinChromLength,
+        chromlen=lambda wildcards: str(getMinChromLength(kmc_log_path=os.path.join(dir.out.kmc, f"{wildcards.sample}", f"{wildcards.sample}_kmcLOG.txt"), sample=wildcards.sample,auto=AUTO)),
     conda:
         os.path.join(dir.env, "scripts.yaml")
     resources:
@@ -104,3 +129,4 @@ rule estimate_sr_coverage:
     threads: config.resources.sml.cpu
     script:
         os.path.join(dir.scripts, "sr_coverage.py")
+
diff --git a/hybracter/workflow/rules/processing/estimate_chromosome.smk b/hybracter/workflow/rules/processing/estimate_chromosome.smk
new file mode 100644
index 00000000..b97cfbb6
--- /dev/null
+++ b/hybracter/workflow/rules/processing/estimate_chromosome.smk
@@ -0,0 +1,78 @@
+"""
+@richardstoeckl with a nice implementation
+https://github.com/gbouras13/hybracter/issues/90
+"""
+rule sanitiseInputReads:
+    input:
+        fastq=get_input_lr_fastqs
+    output:
+        fastq=temp(os.path.join(dir.out.qc, "{sample}_sanitised.fastq.gz"))
+    conda:
+        os.path.join(dir.env, "seqkit.yaml")
+    threads:
+        1
+    shell:
+        """
+        seqkit replace -p "\t" -r '_' {input.fastq} -o {output.fastq}
+        """
+
+rule kmc:
+    """
+    need params.prefix so kmc writes into the output dir
+    """
+    input:
+        fastq = rules.sanitiseInputReads.output.fastq,
+    output:
+        kmcLOG = os.path.join(dir.out.kmc,"{sample}", "{sample}_kmcLOG.txt"),
+        dir = directory(os.path.join(dir.out.kmc,"{sample}")),
+    conda:
+        os.path.join(dir.env, "kmc.yaml")
+    threads:
+        config.resources.med.cpu
+    resources:
+        mem_mb=config.resources.med.mem,
+        mem=str(config.resources.med.mem) + "MB",
+        time=config.resources.med.time,
+    params:
+        prefix=os.path.join(dir.out.kmc,"{sample}","{sample}")
+    shell:
+        "kmc -k21 -ci10 -t{threads} -v -fq {input.fastq} {params.prefix} {output.dir} > {output.kmcLOG} 2>&1"
+
+rule write_chrom_size:
+    """
+    need params.prefix so kmc writes into the output dir
+    """
+    input:
+        kmcLOG = os.path.join(dir.out.kmc,"{sample}", "{sample}_kmcLOG.txt")
+    output:
+        chrom_size = os.path.join(dir.out.chrom_size,"{sample}", "{sample}_kmc_estimated_chrom_size.txt"),
+    conda:
+        os.path.join(dir.env, "scripts.yaml")
+    threads:
+        config.resources.sml.cpu
+    resources:
+        mem_mb=config.resources.sml.mem,
+        mem=str(config.resources.sml.mem) + "MB",
+        time=config.resources.sml.time,
+    script:
+        os.path.join(dir.scripts, "write_chrom_size.py")
+
+rule aggr_kmc:
+    """
+    aggregates over all samples
+    """
+    input:
+        expand(os.path.join(dir.out.chrom_size,"{sample}", "{sample}_kmc_estimated_chrom_size.txt"), sample=SAMPLES),
+        expand(os.path.join(dir.out.kmc,"{sample}", "{sample}_kmcLOG.txt"), sample=SAMPLES)
+    output:
+        flag=os.path.join(dir.out.flags, "aggr_kmc.flag"),
+    resources:
+        mem_mb=config.resources.sml.mem,
+        mem=str(config.resources.sml.mem) + "MB",
+        time=config.resources.sml.time,
+    threads: config.resources.sml.cpu
+    shell:
+        """
+        touch {output.flag}
+        """
+
diff --git a/hybracter/workflow/rules/processing/extract_fastas.smk b/hybracter/workflow/rules/processing/extract_fastas.smk
index 35925e5f..e5d224ad 100644
--- a/hybracter/workflow/rules/processing/extract_fastas.smk
+++ b/hybracter/workflow/rules/processing/extract_fastas.smk
@@ -6,10 +6,11 @@ checkpoint check_completeness:
     input:
         fasta=os.path.join(dir.out.assemblies, "{sample}", "assembly.fasta"),
         info=os.path.join(dir.out.assemblies, "{sample}", "assembly_info.txt"),
+        kmc=os.path.join(dir.out.kmc,"{sample}", "{sample}_kmcLOG.txt"),
     output:
         completeness_check=os.path.join(dir.out.completeness, "{sample}.txt"),
     params:
-        min_chrom_length=getMinChromLength,
+        min_chrom_length=lambda wildcards: str(getMinChromLength(kmc_log_path=os.path.join(dir.out.kmc, f"{wildcards.sample}", f"{wildcards.sample}_kmcLOG.txt"), sample=wildcards.sample,auto=AUTO)),
     conda:
         os.path.join(dir.env, "scripts.yaml")
     resources:
@@ -30,11 +31,12 @@ rule extract_chromosome_complete:
         fasta=os.path.join(dir.out.assemblies, "{sample}", "assembly.fasta"),
         info=os.path.join(dir.out.assemblies, "{sample}", "assembly_info.txt"),
         completeness_check=os.path.join(dir.out.completeness, "{sample}.txt"),
+        kmc=os.path.join(dir.out.kmc,"{sample}", "{sample}_kmcLOG.txt")
     output:
         fasta=os.path.join(dir.out.chrom_pre_polish, "{sample}_chromosome.fasta"),
         ignore_list=os.path.join(dir.out.chrom_pre_polish, "{sample}_ignore_list.txt"),
     params:
-        min_chrom_length=getMinChromLength,
+        min_chrom_length=lambda wildcards: str(getMinChromLength(kmc_log_path=os.path.join(dir.out.kmc, f"{wildcards.sample}", f"{wildcards.sample}_kmcLOG.txt"), sample=wildcards.sample,auto=AUTO)),
         polypolish_flag=False,
     conda:
         os.path.join(dir.env, "scripts.yaml")
@@ -59,7 +61,7 @@ rule copy_flye_intermediate_chrom_assembly:
             dir.out.intermediate_assemblies, "{sample}", "{sample}_flye.fasta"
         ),
     params:
-        min_chrom_length=getMinChromLength,
+        min_chrom_length=lambda wildcards: str(getMinChromLength(kmc_log_path=os.path.join(dir.out.kmc, f"{wildcards.sample}", f"{wildcards.sample}_kmcLOG.txt"), sample=wildcards.sample,auto=AUTO)),
     conda:
         os.path.join(dir.env, "scripts.yaml")
     resources:
@@ -104,10 +106,11 @@ rule extract_incomplete:
     input:
         fasta=os.path.join(dir.out.assemblies, "{sample}", "assembly.fasta"),
         completeness_check=os.path.join(dir.out.completeness, "{sample}.txt"),
+        kmc=os.path.join(dir.out.kmc,"{sample}", "{sample}_kmcLOG.txt")
     output:
         fasta=os.path.join(dir.out.incomp_pre_polish, "{sample}.fasta"),
     params:
-        min_chrom_length=getMinChromLength,
+        min_chrom_length=lambda wildcards: str(getMinChromLength(kmc_log_path=os.path.join(dir.out.kmc, f"{wildcards.sample}", f"{wildcards.sample}_kmcLOG.txt"), sample=wildcards.sample,auto=AUTO)),
     conda:
         os.path.join(dir.env, "scripts.yaml")
     resources:
diff --git a/hybracter/workflow/rules/processing/qc.smk b/hybracter/workflow/rules/processing/qc.smk
index 114d6396..d91c94d3 100644
--- a/hybracter/workflow/rules/processing/qc.smk
+++ b/hybracter/workflow/rules/processing/qc.smk
@@ -5,6 +5,7 @@ rule filtlong:
     """
     input:
         fastq=get_input_lr_fastqs,
+        kmc=os.path.join(dir.out.kmc,"{sample}", "{sample}_kmcLOG.txt")
     output:
         fastq=temp(os.path.join(dir.out.qc, "{sample}_filt.fastq.gz")),
         version=os.path.join(dir.out.versions, "{sample}", "filtlong.version"),
@@ -18,7 +19,7 @@ rule filtlong:
     params:
         qual=config.args.min_quality,
         length=config.args.min_length,
-        target_bases=getTargetBases,
+        target_bases=lambda wildcards: str(getTargetBases(kmc_log_path=os.path.join(dir.out.kmc,f"{wildcards.sample}", f"{wildcards.sample}_kmcLOG.txt"),sample=wildcards.sample, auto=AUTO, target_depth=SUBSAMPLE_DEPTH)),
     benchmark:
         os.path.join(dir.out.bench, "filtlong", "{sample}.txt")
     log:
diff --git a/hybracter/workflow/rules/processing/remove_contaminants_qc.smk b/hybracter/workflow/rules/processing/remove_contaminants_qc.smk
index e1ac23cf..dfee529f 100644
--- a/hybracter/workflow/rules/processing/remove_contaminants_qc.smk
+++ b/hybracter/workflow/rules/processing/remove_contaminants_qc.smk
@@ -32,7 +32,7 @@ rule host_removal_mapping_single:
     """Map reads to host and return unmapped reads"""
     input:
         index=os.path.join(dir.out.contaminant_index, "host.index"),
-        fastq=get_input_lr_fastqs,
+        fastq=get_input_lr_fastqs(sample=wildcards.sample),
     output:
         r1=os.path.join(
             dir.out.contaminant_removal, "{sample}", "{sample}.host_rm.fastq.gz"
@@ -98,6 +98,7 @@ rule filtlong:
         fastq=os.path.join(
             dir.out.contaminant_removal, "{sample}", "{sample}.host_rm.fastq.gz"
         ),
+        kmc=os.path.join(dir.out.kmc,"{sample}", "{sample}_kmcLOG.txt")
     output:
         fastq=temp(os.path.join(dir.out.qc, "{sample}_filt.fastq.gz")),
         version=os.path.join(dir.out.versions, "{sample}", "filtlong.version"),
@@ -111,13 +112,14 @@ rule filtlong:
     params:
         qual=config.args.min_quality,
         length=config.args.min_length,
+        target_bases=lambda wildcards: str(getTargetBases(kmc_log_path=os.path.join(dir.out.kmc,f"{wildcards.sample}", f"{wildcards.sample}_kmcLOG.txt"),sample=wildcards.sample, auto=AUTO, target_depth=SUBSAMPLE_DEPTH)),
     benchmark:
         os.path.join(dir.out.bench, "filtlong", "{sample}.txt")
     log:
         os.path.join(dir.out.stderr, "filtlong", "{sample}.log"),
     shell:
         """
-        filtlong --min_mean_q {params.qual} --min_length {params.length} {input.fastq} | pigz > {output.fastq} 2> {log}
+        filtlong --target_bases {params.target_bases}  --min_mean_q {params.qual} --min_length {params.length} {input.fastq} | pigz > {output.fastq} 2> {log}
         filtlong --version > {output.version}
         """
 
diff --git a/hybracter/workflow/rules/reorientation/dnaapler.smk b/hybracter/workflow/rules/reorientation/dnaapler.smk
index 0ee5a148..1faa1831 100644
--- a/hybracter/workflow/rules/reorientation/dnaapler.smk
+++ b/hybracter/workflow/rules/reorientation/dnaapler.smk
@@ -31,6 +31,6 @@ rule dnaapler:
         os.path.join(dir.out.stderr, "dnaapler", "{sample}.log"),
     shell:
         """
-        dnaapler all -i {input.fasta} -o {params.dir} --ignore {input.ignore_list} -p {wildcards.sample} -t {threads} -a nearest --db dnaa,repa -f 2> {log}
+        dnaapler all -i {input.fasta} -o {params.dir} --ignore {input.ignore_list} -p {wildcards.sample} -t {threads} -a nearest --db dnaa,repa,cog1474 -f 2> {log}
         dnaapler --version > {output.version}
         """
diff --git a/hybracter/workflow/rules/reorientation/dnaapler_no_medaka.smk b/hybracter/workflow/rules/reorientation/dnaapler_no_medaka.smk
index f90b209b..502b5e36 100644
--- a/hybracter/workflow/rules/reorientation/dnaapler_no_medaka.smk
+++ b/hybracter/workflow/rules/reorientation/dnaapler_no_medaka.smk
@@ -28,7 +28,7 @@ rule dnaapler_no_medaka:
         os.path.join(dir.out.stderr, "dnaapler", "{sample}.log"),
     shell:
         """
-        dnaapler all -i {input.fasta} -o {params.dir} --ignore {input.ignore_list} -p {wildcards.sample} -t {threads} -a nearest --db dnaa,repa -f 2> {log}
+        dnaapler all -i {input.fasta} -o {params.dir} --ignore {input.ignore_list} -p {wildcards.sample} -t {threads} -a nearest --db dnaa,repa,cog1474 -f 2> {log}
         dnaapler --version > {output.version}
         """
 
diff --git a/hybracter/workflow/scripts/lr_coverage.py b/hybracter/workflow/scripts/lr_coverage.py
new file mode 100644
index 00000000..1a37dbc9
--- /dev/null
+++ b/hybracter/workflow/scripts/lr_coverage.py
@@ -0,0 +1,47 @@
+#!/usr/bin/env python3
+
+
+import pandas as pd
+import sys
+
+
+def get_sum_len(file_path):
+    """
+    get sum_len from seqkit output
+    """
+
+    # Define the column names based on the file format
+
+    # Read the file into a pandas DataFrame
+    df = pd.read_csv(file_path, sep="\t")
+    sum_len_value = df["sum_len"].values[0]
+
+    return sum_len_value
+
+
+def get_quick_coverage_estimate_long(
+    long_seqkit, min_chrom_length, min_bases, coverage_output_path
+):
+    """
+    takes seqkit input
+    """
+
+    long_bases = get_sum_len(long_seqkit)
+
+    # number of bases is over the min_bases
+    if long_bases > int(min_bases):
+        coverage = round(long_bases / int(min_chrom_length))
+
+        # write to file
+        with open(coverage_output_path, "w") as file:
+            file.write(str(coverage))
+    else: # number of bases is less than min bases
+        sys.exit(f"The total number of bases of long-read sequencing is {long_bases}, which is lower than the required minimum value {min_bases} based on the chromosome size multiplied by the --min_depth parameter.\n Please check --min_depth and/or get more long-read sequencing data! Hybracter will stop here.")
+
+
+get_quick_coverage_estimate_long(
+    snakemake.input.long_bases,
+    snakemake.params.min_chrom_length,
+    snakemake.params.min_bases,
+    snakemake.output.lr_coverage,
+)
diff --git a/hybracter/workflow/scripts/write_chrom_size.py b/hybracter/workflow/scripts/write_chrom_size.py
new file mode 100644
index 00000000..bfd999aa
--- /dev/null
+++ b/hybracter/workflow/scripts/write_chrom_size.py
@@ -0,0 +1,24 @@
+#!/usr/bin/env python3
+
+import re
+
+
+def write_chrom_size(
+    kmc_log_path, output_file
+):
+
+    with open(kmc_log_path, "r") as file:
+        for line in file:
+            if "No. of unique counted k-mers" in line:
+                # keep 80% of kmers as lower bound for chromosome
+                chrom_size = int(float(re.search(r"No. of unique counted k-mers\s*:\s*([\d\.eE+-]+)", line).group(1)) * 0.8)
+
+    # write to file
+    with open(output_file, "w") as file:
+        file.write(str(chrom_size))
+
+
+write_chrom_size(
+    snakemake.input.kmcLOG,
+    snakemake.output.chrom_size,
+)


Test patch:
diff --git a/.github/workflows/testing.yml b/.github/workflows/testing.yml
index 85e76ad2..0e361fcf 100644
--- a/.github/workflows/testing.yml
+++ b/.github/workflows/testing.yml
@@ -13,7 +13,7 @@ jobs:
 
       strategy:
         matrix:
-          os: [macos-latest]
+          os: [macos-13] #macos-latest/macos-14 is M1 - porechop-abi etc don't have M1 package builds, so fails CI
           python-version: ["3.11"]
 
       steps:
@@ -22,20 +22,19 @@ jobs:
             fetch-depth: 0
 
         # Setup env
-        - uses: conda-incubator/setup-miniconda@v2
+        - uses: conda-incubator/setup-miniconda@v3
           with:
             activate-environment: hybracterENV
             environment-file: build/environment.yaml
             python-version: ${{ matrix.python-version }}
             auto-activate-base: false
-            miniforge-variant: Mambaforge
-            channels: conda-forge,bioconda,defaults,anaconda
+            channel-priority: strict
             auto-update-conda: true
 
         - name: "Setup hybracter on ${{ matrix.os }} for Python ${{ matrix.python-version }}"
           shell: bash -l {0}
           run: |
-            mamba install python=${{ matrix.python-version }}
+            conda install python=${{ matrix.python-version }}
             python -m pip install --upgrade pip
             pip install .
         
@@ -66,20 +65,20 @@ jobs:
           fetch-depth: 0
 
       # Setup env
-      - uses: conda-incubator/setup-miniconda@v2
+      - uses: conda-incubator/setup-miniconda@v3
         with:
           activate-environment: hybracterENV
           environment-file: build/environment.yaml
           python-version: ${{ matrix.python-version }}
           auto-activate-base: false
-          miniforge-variant: Mambaforge
+          channel-priority: strict
           channels: conda-forge,bioconda,defaults,anaconda
           auto-update-conda: true
 
       - name: "Setup hybracter on ${{ matrix.os }} for Python ${{ matrix.python-version }}"
         shell: bash -l {0}
         run: |
-          mamba install python=${{ matrix.python-version }}
+          conda install python=${{ matrix.python-version }}
           python -m pip install --upgrade pip
           pip install .
       
diff --git a/hybracter/test_data/Fastqs/test_long_reads_subsampled.fastq.gz b/hybracter/test_data/Fastqs/test_long_reads_subsampled.fastq.gz
new file mode 100644
index 00000000..1c0a2477
Binary files /dev/null and b/hybracter/test_data/Fastqs/test_long_reads_subsampled.fastq.gz differ
diff --git a/hybracter/test_data/test_hybrid_input_auto.csv b/hybracter/test_data/test_hybrid_input_auto.csv
new file mode 100644
index 00000000..5504b230
--- /dev/null
+++ b/hybracter/test_data/test_hybrid_input_auto.csv
@@ -0,0 +1,2 @@
+Sample1,hybracter/test_data/Fastqs/test_long_reads.fastq.gz,hybracter/test_data/Fastqs/test_short_reads_R1.fastq.gz,hybracter/test_data/Fastqs/test_short_reads_R2.fastq.gz
+Sample2,hybracter/test_data/Fastqs/test_long_reads_subsampled.fastq.gz,hybracter/test_data/Fastqs/test_short_reads_R1.fastq.gz,hybracter/test_data/Fastqs/test_short_reads_R2.fastq.gz
diff --git a/hybracter/test_data/test_long_input_auto.csv b/hybracter/test_data/test_long_input_auto.csv
new file mode 100644
index 00000000..b7dc2bce
--- /dev/null
+++ b/hybracter/test_data/test_long_input_auto.csv
@@ -0,0 +1,2 @@
+Sample1,hybracter/test_data/Fastqs/test_long_reads.fastq.gz
+Sample2,hybracter/test_data/Fastqs/test_long_reads_subsampled.fastq.gz
diff --git a/hybracter/workflow/test_hybrid.smk b/hybracter/workflow/test_hybrid.smk
index 352ac77a..edc36055 100644
--- a/hybracter/workflow/test_hybrid.smk
+++ b/hybracter/workflow/test_hybrid.smk
@@ -63,18 +63,24 @@ FLYE_MODEL = config.args.flyeModel
 LOGIC = config.args.logic
 DEPTH_FILTER = config.args.depth_filter
 SUBSAMPLE_DEPTH = config.args.subsample_depth
+MIN_DEPTH = config.args.min_depth
+AUTO = config.args.auto
 MAC = config.args.mac
 
 # MAC medaka
 
-new_models = ['r1041_e82_400bps_sup_v5.0.0',
-    'r1041_e82_400bps_hac_v5.0.0',
-    'r1041_e82_400bps_hac_v4.3.0',
-    'r1041_e82_400bps_sup_v4.3.0']
+new_models = [
+    "r1041_e82_400bps_sup_v5.0.0",
+    "r1041_e82_400bps_hac_v5.0.0",
+    "r1041_e82_400bps_hac_v4.3.0",
+    "r1041_e82_400bps_sup_v4.3.0",
+]
 
 if MAC:
     if MEDAKA_MODEL in new_models:
-        print(f"{MEDAKA_MODEL} is not available in medaka v1.8.0 as it is too new. If you want this model, try Hybracter on a Linux machine.")
+        print(
+            f"{MEDAKA_MODEL} is not available in medaka v1.8.0 as it is too new. If you want this model, try Hybracter on a Linux machine."
+        )
         print(f"Changing the medaka model to r1041_e82_400bps_sup_v4.2.0.")
         MEDAKA_MODEL = "r1041_e82_400bps_sup_v4.2.0"
 
@@ -83,15 +89,19 @@ if MAC:
 dictReads = {}
 dictReads["Sample1"] = {}
 dictReads["Sample1"]["LR"] = os.path.join(dir.test_fastqs, "test_long_reads.fastq.gz")
-dictReads["Sample1"]["MinChromLength"] = 50000
+if not AUTO:
+    dictReads["Sample1"]["MinChromLength"] = 50000
 dictReads["Sample1"]["TargetBases"] = SUBSAMPLE_DEPTH * 50000
+dictReads["Sample1"]["MinBases"] = MIN_DEPTH * 50000
 dictReads["Sample1"]["R1"] = os.path.join(dir.test_fastqs, "test_short_reads_R1.fastq.gz")
 dictReads["Sample1"]["R2"] = os.path.join(dir.test_fastqs, "test_short_reads_R2.fastq.gz")
 
 dictReads["Sample2"] = {}
 dictReads["Sample2"]["LR"] = os.path.join(dir.test_fastqs, "test_long_reads.fastq.gz")
-dictReads["Sample2"]["MinChromLength"] = 100000
+if not AUTO:
+    dictReads["Sample2"]["MinChromLength"] = 100000
 dictReads["Sample2"]["TargetBases"] = SUBSAMPLE_DEPTH * 100000
+dictReads["Sample2"]["MinBases"] = MIN_DEPTH * 100000
 dictReads["Sample2"]["R1"] = os.path.join(dir.test_fastqs, "test_short_reads_R1.fastq.gz")
 dictReads["Sample2"]["R2"] = os.path.join(dir.test_fastqs, "test_short_reads_R2.fastq.gz")
 
@@ -103,6 +113,9 @@ SAMPLES = ["Sample1", "Sample2"]
 # Import rules and functions
 ##############################
 
+# kmc - needs to be included due to the chckpointing
+include: os.path.join("rules", "processing", "estimate_chromosome.smk")
+
 # qc and host
 # depends on whehter --contaminants has been specified and --skip_qc flag activiated
 if config.args.contaminants != "none":  # where --contaminants specified
diff --git a/hybracter/workflow/test_long.smk b/hybracter/workflow/test_long.smk
index 6b47fd1a..f26e94df 100644
--- a/hybracter/workflow/test_long.smk
+++ b/hybracter/workflow/test_long.smk
@@ -63,18 +63,24 @@ FLYE_MODEL = config.args.flyeModel
 LOGIC = config.args.logic
 DEPTH_FILTER = config.args.depth_filter
 SUBSAMPLE_DEPTH = config.args.subsample_depth
+MIN_DEPTH = config.args.min_depth
+AUTO = config.args.auto
 MAC = config.args.mac
 
 # MAC medaka
 
-new_models = ['r1041_e82_400bps_sup_v5.0.0',
-    'r1041_e82_400bps_hac_v5.0.0',
-    'r1041_e82_400bps_hac_v4.3.0',
-    'r1041_e82_400bps_sup_v4.3.0']
+new_models = [
+    "r1041_e82_400bps_sup_v5.0.0",
+    "r1041_e82_400bps_hac_v5.0.0",
+    "r1041_e82_400bps_hac_v4.3.0",
+    "r1041_e82_400bps_sup_v4.3.0",
+]
 
 if MAC:
     if MEDAKA_MODEL in new_models:
-        print(f"{MEDAKA_MODEL} is not available in medaka v1.8.0 as it is too new. If you want this model, try Hybracter on a Linux machine.")
+        print(
+            f"{MEDAKA_MODEL} is not available in medaka v1.8.0 as it is too new. If you want this model, try Hybracter on a Linux machine."
+        )
         print(f"Changing the medaka model to r1041_e82_400bps_sup_v4.2.0.")
         MEDAKA_MODEL = "r1041_e82_400bps_sup_v4.2.0"
 
@@ -89,20 +95,28 @@ if MAC:
 dictReads = {}
 dictReads["Sample1"] = {}
 dictReads["Sample1"]["LR"] = os.path.join(dir.test_fastqs, "test_long_reads.fastq.gz")
-dictReads["Sample1"]["MinChromLength"] = 50000
+if not AUTO:
+    dictReads["Sample1"]["MinChromLength"] = 50000
 dictReads["Sample1"]["TargetBases"] = SUBSAMPLE_DEPTH * 50000
+dictReads["Sample1"]["MinBases"] = MIN_DEPTH * 50000
 
 dictReads["Sample2"] = {}
 dictReads["Sample2"]["LR"] = os.path.join(dir.test_fastqs, "test_long_reads.fastq.gz")
-dictReads["Sample2"]["MinChromLength"] = 100000
+if not AUTO:
+    dictReads["Sample2"]["MinChromLength"] = 100000
 dictReads["Sample2"]["TargetBases"] = SUBSAMPLE_DEPTH * 100000
+dictReads["Sample2"]["MinBases"] = MIN_DEPTH * 100000
 
 SAMPLES = ["Sample1", "Sample2"]
 
+
 ##############################
 # Import rules and functions
 ##############################
 
+# kmc - needs to be included due to the chckpointing
+include: os.path.join("rules", "processing", "estimate_chromosome.smk")
+
 # qc and host
 # depends on whehter --contaminants has been specified and --skip_qc flag activiated
 if config.args.contaminants != "none":  # where --contaminants specified
diff --git a/tests/test_hybracter.py b/tests/test_hybracter.py
index 5e7e0ead..767231dd 100644
--- a/tests/test_hybracter.py
+++ b/tests/test_hybracter.py
@@ -82,7 +82,8 @@ def test_hybracter_t_hybrid_no_medaka(run_mac):
     outdir: Path = "test_hybracter_output"
     cmd = f"hybracter test-hybrid --threads {threads} --output {outdir} --no_medaka"
     if run_mac:
-        cmd += " --mac"
+        cmd = f"hybracter test-hybrid -h"
+        #cmd += " --mac"
     exec_command(cmd)
     remove_directory(outdir)
 
@@ -92,7 +93,8 @@ def test_hybracter_t_subsample_depth(run_mac):
     outdir: Path = "test_hybracter_output"
     cmd = f"hybracter test-hybrid --threads {threads} --output {outdir} --no_medaka --subsample_depth 50"
     if run_mac:
-        cmd += " --mac"
+        cmd = f"hybracter test-hybrid -h"
+        #cmd += " --mac"
     exec_command(cmd)
     remove_directory(outdir)
 
@@ -101,16 +103,18 @@ def test_hybracter_t_depth_filter(run_mac):
     outdir: Path = "test_hybracter_output"
     cmd = f"hybracter test-hybrid --threads {threads} --output {outdir} --no_medaka --depth_filter 2"
     if run_mac:
-        cmd += " --mac"
+        cmd = f"hybracter test-hybrid -h"
+        #cmd += " --mac"
     exec_command(cmd)
     remove_directory(outdir)
 
 def test_hybracter_hybrid_last_no_medaka(run_mac):
-    """hybracter test-hybrid"""
+    """hybracter test-hybrid - don't worry about MAC as no medaka"""
     outdir: Path = "test_hybracter_output"
     cmd = f"hybracter test-hybrid --threads {threads} --output {outdir} --logic last --no_medaka --skip_qc"
     if run_mac:
-        cmd += " --mac"
+        cmd = f"hybracter test-hybrid -h"
+        #cmd += " --mac"
     exec_command(cmd)
     remove_directory(outdir)
 
@@ -164,11 +168,12 @@ def test_hybracter_t_long_depth_filter(run_mac):
     remove_directory(outdir)
 
 def test_hybracter_t_long_no_medaka(run_mac):
-    """hybracter test-long no medaka"""
+    """hybracter test-long no medaka - don't worry about mac"""
     outdir: Path = "test_hybracter_output"
     cmd = f"hybracter test-long --threads {threads} --output {outdir} --no_medaka"
     if run_mac:
-        cmd += " --mac"
+        cmd = f"hybracter test-hybrid -h"
+        #cmd += " --mac"
     exec_command(cmd)
     remove_directory(outdir)
 
@@ -196,6 +201,32 @@ def test_hybracter_hybrid_csv(run_mac):
     exec_command(cmd)
     remove_directory(outdir)
 
+def test_hybracter_hybrid_csv_auto(run_mac):
+    """test hybracter hybrid auto
+
+    Note that in 'auto' mode, it is going to be very hard for incomplete samples to happen
+    This is because the required No. unique counter k-mers to be count is set to 10 
+    Likely to get longer assemblies than 10kmers (in general), done some tests on subsampled data
+    So warn people to be careful using kmc
+
+    Sample 1
+    kmc
+    No. of unique counted k-mers       :        72697
+    which should give a chromosome size of 58157  ( assembly is 70ish)
+
+    Sample 2
+    No. of unique counted k-mers       : 4542  
+    which should give a chromosome size of 3633( assembly is 70ish)
+
+    """
+    outdir: Path = "test_hybracter_output"
+    input_csv: Path = test_data_path / "test_hybrid_input_auto.csv"
+    cmd = f"hybracter hybrid --input {input_csv} --threads {threads} --output {outdir} --databases {db_dir} --auto"
+    if run_mac:
+        cmd += " --mac"
+    exec_command(cmd)
+    remove_directory(outdir)
+
 def test_hybracter_hybrid_csv_skip_qc_not_gzipped(run_mac):
     """test hybracter hybrid with --skip_qc and a long read fastq that isn't gzipped """
     outdir: Path = "test_hybracter_output"
@@ -232,7 +263,8 @@ def test_hybracter_no_pypolca_no_medaka(run_mac):
     input_csv: Path = test_data_path / "test_hybrid_input.csv"
     cmd = f"hybracter hybrid --input {input_csv} --threads {threads} --output {outdir} --no_pypolca --databases {db_dir} --no_medaka "
     if run_mac:
-        cmd += " --mac"
+        cmd = f"hybracter test-hybrid -h"
+        #cmd += " --mac"
     exec_command(cmd)
     remove_directory(outdir)
 
@@ -243,7 +275,8 @@ def test_hybracter_no_medaka(run_mac):
     input_csv: Path = test_data_path / "test_hybrid_input.csv"
     cmd = f"hybracter hybrid --input {input_csv} --threads {threads} --output {outdir} --databases {db_dir} --no_medaka"
     if run_mac:
-        cmd += " --mac"
+        cmd = f"hybracter test-hybrid -h"
+        #cmd += " --mac"
     exec_command(cmd)
     remove_directory(outdir)
 
@@ -268,9 +301,8 @@ def test_hybracter_subsample_depth(run_mac):
     exec_command(cmd)
     remove_directory(outdir)
 
-
-def test_hybracter_hybrid_datadir(run_mac):
-    """test hybracter hybrid with --datadir"""
+def test_hybracter_hybrid_datadir_1dir(run_mac):
+    """test hybracter hybrid with --datadir - 1 dir"""
     outdir: Path = "test_hybracter_output"
     input_csv: Path = test_data_path / "test_hybrid_input_datadir.csv"
     datadir: Path = "hybracter/test_data/Fastqs/"
@@ -280,6 +312,18 @@ def test_hybracter_hybrid_datadir(run_mac):
     exec_command(cmd)
     remove_directory(outdir)
 
+def test_hybracter_hybrid_datadir_2dirs(run_mac):
+    """test hybracter hybrid with --datadir - 2 dirs"""
+    outdir: Path = "test_hybracter_output"
+    input_csv: Path = test_data_path / "test_hybrid_input_datadir.csv"
+    # 2 dirs
+    datadir: Path = "hybracter/test_data/Fastqs/,hybracter/test_data/Fastqs/"
+    cmd = f"hybracter hybrid --input {input_csv} --threads {threads} --output {outdir} --databases {db_dir} --datadir {datadir}"
+    if run_mac:
+        cmd += " --mac"
+    exec_command(cmd)
+    remove_directory(outdir)
+
 """
 long
 """
@@ -294,6 +338,33 @@ def test_hybracter_long(run_mac):
     exec_command(cmd)
     remove_directory(outdir)
 
+def test_hybracter_long_auto(run_mac):
+    """
+
+    Note that in 'auto' mode, it is going to be very hard for incomplete samples to happen
+    This is because the required No. unique counter k-mers to be count is set to 10 
+    Likely to get longer assemblies than 10kmers (in general), done some tests on subsampled data
+    If the chrom size is really low --> suggests that you need more sequencing data :)
+    So warn people to be careful using kmc
+
+    Sample 1
+    kmc
+    No. of unique counted k-mers       :        72697
+    which should give a chromosome size of 58157  ( assembly is 70ish)
+
+    Sample 2
+    No. of unique counted k-mers       : 4542  
+    which should give a chromosome size of 3633( assembly is 70ish)
+
+    """
+    outdir: Path = "test_hybracter_output"
+    input_csv: Path = test_data_path / "test_long_input_auto.csv"
+    cmd = f"hybracter long --input {input_csv} --threads {threads} --output {outdir} --databases {db_dir} --auto"
+    if run_mac:
+        cmd += " --mac"
+    exec_command(cmd)
+    remove_directory(outdir)
+
 def test_hybracter_long_datadir(run_mac):
     """test hybracter long with --datadir"""
     outdir: Path = "test_hybracter_output"
@@ -304,6 +375,7 @@ def test_hybracter_long_datadir(run_mac):
         cmd += " --mac"
     exec_command(cmd)
     remove_directory(outdir)
+    
 
 def test_hybracter_long_no_medaka(run_mac):
     """test hybracter long"""
@@ -311,7 +383,8 @@ def test_hybracter_long_no_medaka(run_mac):
     input_csv: Path = test_data_path / "test_long_input.csv"
     cmd = f"hybracter long --input {input_csv} --threads {threads} --output {outdir} --databases {db_dir} --no_medaka"
     if run_mac:
-        cmd += " --mac"
+        cmd = f"hybracter test-hybrid -h"
+        #cmd += " --mac"
     exec_command(cmd)
     remove_directory(outdir)
 
@@ -342,7 +415,8 @@ def test_hybracter_hybrid_single_no_medaka(run_mac):
     reads2: Path = test_data_path / "Fastqs/test_short_reads_R2.fastq.gz"
     cmd = f"hybracter hybrid-single -l {longreads} -1 {reads1} -2 {reads2} -c 50000 -s Sample1 --threads {threads} --output {outdir} --databases {db_dir} --no_medaka"
     if run_mac:
-        cmd += " --mac"
+        cmd = f"hybracter test-hybrid -h"
+        #cmd += " --mac"
     exec_command(cmd)
     remove_directory(outdir)
 
@@ -364,7 +438,8 @@ def test_hybracter_long_single_no_medaka_subsample_depth(run_mac):
     longreads: Path = test_data_path / "Fastqs/test_long_reads.fastq.gz"
     cmd = f"hybracter long-single -l {longreads} -c 50000 -s Sample1  --threads {threads} --output {outdir} --databases {db_dir} --no_medaka --subsample_depth 50"
     if run_mac:
-        cmd += " --mac"
+        cmd = f"hybracter test-hybrid -h"
+        #cmd += " --mac"
     exec_command(cmd)
     remove_directory(outdir)
 
@@ -400,17 +475,42 @@ def test_version():
     exec_command(cmd)
 
 
-# class errors(unittest.TestCase):
-#     def test_no_db(self, run_mac):
-#         """test hybracter with no db"""
-#         with self.assertRaises(RuntimeError):
-#             outdir: Path = "test_hybracter_output"
-#             empty_db: Path = "tests/empty_db"
-#             cmd = f"hybracter test-hybrid --threads {threads} --output {outdir} --database {empty_db} --no_pypolca"
-#             if run_mac:
-#                 cmd += " --mac"
-#             exec_command(cmd)
-#             remove_directory(outdir)
+
+class TestErrors:
+    @pytest.fixture(autouse=True)
+    def setup(self):
+        self.threads = threads  
+
+    def test_no_db(self, run_mac):
+        """test hybracter with no db"""
+        with pytest.raises(RuntimeError):
+            outdir: Path = "test_hybracter_output"
+            empty_db: Path = "tests/empty_db"
+            cmd = f"hybracter test-hybrid --threads {threads} --output {outdir} --database {empty_db} --no_pypolca"
+            if run_mac:
+                cmd += " --mac"
+            exec_command(cmd)
+            remove_directory(outdir)
+
+    def test_hybracter_t_hybrid_md(run_mac):
+        """hybracter test-hybrid with --min_depth 1000"""
+        with pytest.raises(RuntimeError):
+            outdir: Path = "test_hybracter_output"
+            cmd = f"hybracter test-hybrid --threads {threads} --output {outdir} --min_depth 1000"
+            if run_mac:
+                cmd += " --mac"
+            exec_command(cmd)
+            remove_directory(outdir)
+
+    def test_hybracter_t_long_md(run_mac):
+        """hybracter test-long --min_depth 1000 also test """
+        with pytest.raises(RuntimeError):
+            outdir: Path = "test_hybracter_output"
+            cmd = f"hybracter test-long --threads {threads} --output {outdir} --min_depth 1000"
+            if run_mac:
+                cmd += " --mac"
+            exec_command(cmd)
+            remove_directory(outdir)
 
 
 # cleanup