shell:
    """
     [ -d logs/ ] || mkdir -p logs/
     cd logs/
     for RULE in {RULES}; do
      [ -d ${{RULE}}/ ] || mkdir -p ${{RULE}}/
     done
    """

run:
    for old,new in zip(metadata.Old_name,metadata.New_name):
        shell("[ -f {new} ] || mv {old} {new}".format(old=old, new=new))

run:
    for old,new in zip(metadata.Old_name,metadata.New_name):
        shell("mv {old} {new}".format(old=old, new=new))

shell:
    """
    set +u
    {params.conda_activate}
    {params.PERL5LIB}
    set -u

      {params.program} --outdir {params.out_dir}/ \
         --threads {params.threads} {input.forward} {input.rev}

    """

shell:
    """
    set +u
    {params.conda_activate}
    {params.PERL5LIB}
    set -u

      {params.program} \
          --interactive \
          -f {params.out_dir} \
          -o {params.out_dir}
    """

shell:
    """
    set +u
    {params.conda_activate}
    set -u

    qiime tools import \
         --type 'SampleData[SequencesWithQuality]' \
         --input-path {input.manifest_file} \
         --output-path {output} \
         --input-format SingleEndFastqManifestPhred33
    """

shell:
    """
    set +u
    {params.conda_activate}
    set -u

    qiime cutadapt trim-single \
         --i-demultiplexed-sequences {input} \
         --p-cores {params.cores} \
         --p-front {params.forward_primer} \
         --o-trimmed-sequences {output} \
         --verbose

    """

shell:
    """
    {params.conda_activate}

     [ -d {params.out_dir} ] ||  mkdir -p {params.out_dir}
     # Merge reads then delete unnecessary files
     {params.program} \
        -f {input.forward} \
        -r {input.rev} \
        -j {params.threads} \
        -o {params.out_dir}/{wildcards.sample} \
        -m {params.max} \
        -n {params.min} \
        -t {params.min_trim} > {log} 2>&1


     rm -rf \
       {params.out_dir}/{wildcards.sample}.discarded.fastq \
       {params.out_dir}/{wildcards.sample}.unassembled.forward.fastq \
       {params.out_dir}/{wildcards.sample}.unassembled.reverse.fastq 

     mv {params.out_dir}/{wildcards.sample}.assembled.fastq {params.out_dir}/{wildcards.sample}.fastq

     # gzip to save memory

     gzip {params.out_dir}/{wildcards.sample}.fastq

   """

shell:
    """
    set +u
    {params.conda_activate}
    set -u

    qiime tools import \
         --type 'SampleData[SequencesWithQuality]' \
         --input-path {input.manifest_file} \
         --output-path {output} \
         --input-format SingleEndFastqManifestPhred33
    """

shell:
    """
    set +u
    {params.conda_activate}
    set -u

    qiime cutadapt trim-single \
         --i-demultiplexed-sequences {input} \
         --p-cores {params.cores} \
         --p-front {params.forward_primer} \
         --o-trimmed-sequences {output} \
         --verbose
    """

shell:
    """
    set +u
    {params.conda_activate}
    set -u

    qiime tools import \
         --type 'SampleData[PairedEndSequencesWithQuality]' \
         --input-path {input.manifest_file} \
         --output-path {output} \
         --input-format PairedEndFastqManifestPhred33
    """

shell:
    """
    set +u
    {params.conda_activate}
    set -u

    qiime cutadapt trim-single \
         --i-demultiplexed-sequences {input} \
         --p-cores {params.cores} \
         --p-front {params.forward_primer} \
         --o-trimmed-sequences {output} \
         --verbose
    """

shell:
    """
    set +u
    {params.conda_activate}
    set -u

    qiime demux summarize \
        --p-n 10000 \
        --i-data {input} \
        --o-visualization {output}
"""

shell:
    """
    set +u
    {params.conda_activate}
    set -u

    qiime demux summarize \
        --p-n 10000 \
        --i-data {input} \
        --o-visualization {output}
    """

shell:
    """
    set +u
    {params.conda_activate}
    set -u

    qiime demux summarize \
        --p-n 10000 \
        --i-data {input} \
        --o-visualization {output}
    """

shell:
    """
    set +u
    {params.conda_activate}
    set -u

    MODE={params.mode}

    if [ ${{MODE}} == "paired" ];then

        # Paired end
        qiime dada2 denoise-paired \
            --i-demultiplexed-seqs {input} \
            --o-table {output.table} \
            --o-representative-sequences {output.rep_seqs} \
            --o-denoising-stats {output.stats} \
            --p-trunc-len-f  {params.trun_len_forward} \
            --p-trunc-len-r {params.trun_len_reverse} \
            --p-trim-left-f {params.trim_len_forward} \
            --p-trim-left-r {params.trim_len_reverse} \
            --p-max-ee-f {params.max_forward_err} \
            --p-max-ee-r {params.max_reverse_err} \
            --p-n-threads {params.threads} 

    else

        # Single end
        qiime dada2 denoise-single \
            --i-demultiplexed-seqs {input} \
            --o-table {output.table} \
            --o-representative-sequences {output.rep_seqs} \
            --o-denoising-stats {output.stats} \
            --p-trunc-len  {params.trun_len_forward} \
            --p-trim-left {params.trim_len_forward} \
            --p-max-ee {params.max_forward_err} \
            --p-n-threads {params.threads}

    fi

    """

        shell:
            """
            set +u
            {params.conda_activate}
            set -u

            # Initial quality filtering process based on quality scores
            qiime quality-filter q-score \
              --i-demux {input} \
              --o-filtered-sequences {output.filtered_reads} \
              --o-filter-stats {output.filter_stats}

            # Tabulate the filter statistics
            qiime metadata tabulate \
	           --m-input-file {output.filter_stats} \
 	           --o-visualization {output.filter_stats_viz}


            # # Next, the Deblur workflow is applied using the qiime deblur denoise-16S method.
            # This method requires one parameter that is used in quality filtering,
            # --p-trim-length n which truncates the sequences at position n.
            #  In general, the Deblur developers recommend setting this value to a length 
            # where the median quality score begins to drop too low

            qiime deblur denoise-16S \
              --i-demultiplexed-seqs {output.filtered_reads} \
              --p-trim-length {params.trunc_length} \
              --o-representative-sequences {params.rep_seqs} \
              --o-table {ouput.table} \
              --p-sample-stats \
              --o-stats {output.stats}
            """

    shell:
        """
        set +u
        {params.conda_activate}
        set -u

        qiime feature-table summarize \
	       --i-table {input} \
	       --o-visualization {output}
        """

shell:
    """
    set +u
    {params.conda_activate}
    set -u

    qiime feature-table tabulate-seqs \
       --i-data {input} \
       --o-visualization {output}
    """

shell:
    """
    set +u
    {params.conda_activate}
    set -u

    # Visualize dada2 denoise stats
    qiime metadata tabulate \
      --m-input-file {input} \
      --o-visualization {output}
"""

shell:
    """
    set +u
    {params.conda_activate}
    set -u

     # Visualize deblur stats
     qiime deblur visualize-stats \
         --i-deblur-stats {output.stats} \
         --o-visualization {output}
    """

shell:
    """
    set +u
    {params.conda_activate}
    set -u

     # Assign taxonomy
     qiime feature-classifier classify-sklearn \
       --i-classifier {input.classifier} \
       --i-reads {input.rep_seqs} \
       --o-classification {output.raw} \
       --p-n-jobs {params.threads}

     # Tabulate taxonomy

     qiime metadata tabulate \
       --m-input-file {output.raw} \
       --o-visualization {output.viz}
    """

shell:
    """
    set +u
    {params.conda_activate}
    set -u

    # Remove singletons
    qiime feature-table filter-features \
      --i-table {input} \
      --p-min-frequency 2 \
      --o-filtered-table {output.table_raw}

    qiime feature-table summarize \
      --i-table {output.table_raw} \
      --o-visualization {output.table_viz}
    """

shell:
    """
    set +u
    {params.conda_activate}
    set -u

    # Filter out non-target assigments
    if [ {params.amplicon} == "ITS" ]; then

        # Retain only Fungi sequences
        qiime taxa filter-table \
          --i-table {input.table} \
          --i-taxonomy  {input.taxonomy} \
          --p-include  {params.taxa2exclude} \
          --o-filtered-table {output.table_raw}


    else

        # Filter out non-target assigments
        qiime taxa filter-table \
          --i-table {input.table} \
          --i-taxonomy  {input.taxonomy} \
          --p-exclude  {params.taxa2exclude} \
          --o-filtered-table {output.table_raw}

    fi

    # To figure out the total number of sequences ("Total freqency") 
    # to be used to determine the minuminum frequency for filtering out
    # rare taxa. to calculate the multiply the total number of sequences
    # by 0.005
    qiime feature-table summarize \
      --i-table {output.table_raw} \
      --o-visualization {output.table_viz}
    """

shell:
    """
    set +u
    {params.conda_activate}
    set -u

    # Removing rare otus / features with abundance less the 0.005%
    qiime feature-table filter-features \
      --i-table {input} \
      --p-min-frequency {params.minumum_frequency} \
      --o-filtered-table {output.table_raw}

    qiime feature-table summarize \
      --i-table {output.table_raw} \
      --o-visualization {output.table_viz}
    """

shell:
    """
    set +u
    {params.conda_activate}
    set -u


   # Samples bar plot
   qiime taxa barplot \
     --i-table {input.table} \
     --i-taxonomy {input.taxonomy} \
     --m-metadata-file {input.metadata} \
     --o-visualization  {output}
   """

shell:
    """
    set +u
    {params.conda_activate}
    set -u

    # Group feature table by group in metadata file
    qiime feature-table group \
        --i-table  {input.table}  \
        --p-axis sample \
        --m-metadata-file {input.metadata} \
        --m-metadata-column '{params.category}' \
        --p-mode {params.mode} \
        --o-grouped-table {output.grouped_table}

    # Grouped bar plot
    qiime taxa barplot \
      --i-table {output.grouped_table} \
      --i-taxonomy {input.taxonomy} \
      --m-metadata-file {params.metadata} \
      --o-visualization  {output.bar_plot}
  """

shell:
    """
    set +u
    {params.conda_activate}
    set -u

    qiime diversity core-metrics-phylogenetic \
       --p-sampling-depth {params.depth} \
       --i-table {input.table} \
       --i-phylogeny {input.tree} \
       --m-metadata-file {input.metadata} \
       --p-n-jobs-or-threads 'auto' \
       --output-dir core_diversity/  && \
       mv core_diversity/* {diversity_dir}/ && \
       rm -rf core_diversity/
    """

shell:
    """
    set +u
    {params.conda_activate}
    set -u

    qiime diversity alpha-rarefaction \
       --p-max-depth {params.depth} \
       --i-table {input.table} \
       --i-phylogeny {input.tree} \
       --m-metadata-file {input.metadata} \
       --o-visualization {output}
    """

shell:
    """
    set +u
    {params.conda_activate}
    set -u

    for metric in {alpha_diversity_metrics}; do

        qiime diversity alpha-group-significance \
           --i-alpha-diversity  {diversity_dir}/${{metric}}_vector.qza \
           --m-metadata-file {input.metadata} \
           --o-visualization {diversity_dir}/alpha_${{metric}}_significance.qzv

    done
    """

shell:
    """
    set +u
    {params.conda_activate}
    set -u

    for distance in {distance_matrices}; do

        qiime diversity beta-group-significance \
           --i-distance-matrix  {diversity_dir}/${{distance}}_distance_matrix.qza \
           --m-metadata-file {input.metadata} \
           --m-metadata-column {params.category} \
           --o-visualization {diversity_dir}/beta_${{distance}}_significance.qzv

    done

    """

shell:
    """
    set +u
    {params.conda_activate}
    set -u

    TAXON_LEVELS=(2 3 4 5 6)

    for TAXON_LEVEL in ${{TAXON_LEVELS[*]}}; do

        # Collapse ASV table at a taxonomy level of interest
        qiime taxa collapse \
            --i-table {input.table} \
            --i-taxonomy {input.taxonomy} \
            --p-level ${{TAXON_LEVEL}} \
            --o-collapsed-table {params.out_dir}/L${{TAXON_LEVEL}}-filtered_table.qza

    done
    """

shell:
    "cp {input} {params.out_dir}/  && "
    "mv {params.out_dir}/{params.basename} {output}"

shell:
    """
    set +u
    {params.conda_activate}
    set -u

    qiime composition add-pseudocount \
        --i-table {input} \
        --o-composition-table {output}

    """

shell:
    """
    set +u
    {params.conda_activate}
    set -u

    # Apply ANCOM to identify ASV/OTUs that differ in abundance
    qiime composition ancom \
        --i-table {input.table} \
        --m-metadata-file {input.metadata} \
        --m-metadata-column {params.category} \
        --o-visualization {output}

    """

shell:
    """
    set +u
    {params.conda_activate}
    set -u

    # Export feature table
    qiime tools export \
       --input-path {input.feature_table} \
       --output-path {params.out_dir}/


    # Export representative sequences
    qiime tools export --input-path {input.rep_seqs} --output-path {params.out_dir}/ && \
    mv {params.out_dir}/dna-sequences.fasta  {output.rep_seqs}

    # Export taxonomy
    qiime tools export \
         --input-path {input.taxonomy_table} \
         --output-path {params.out_dir}/


    # ---------------------- Add taxonomy to feature table ------------------------ #

    # Creating a TSV BIOM table
    biom convert \
        -i {params.out_dir}/feature-table.biom \
        -o {params.out_dir}/feature-table.tsv \
        --to-tsv

    # Next, we’ll need to modify the exported taxonomy file’s header before using it with BIOM software.

    # Before modifying that file, make a copy:
    cp {params.out_dir}/taxonomy.tsv {params.out_dir}/biom-taxonomy.tsv

    # Change the first line of biom-taxonomy.tsv (i.e. the header) to this:
    # Note that you’ll need to use tab characters in the header since this is a TSV file.
    #OTUID	taxonomy	confidence   

    (echo "#OTUID	taxonomy	confidence"; sed -e '1d' {params.out_dir}/biom-taxonomy.tsv) \
     > {params.out_dir}/tmp.tsv && \
     rm -rf {params.out_dir}/biom-taxonomy.tsv && \
     mv {params.out_dir}/tmp.tsv {params.out_dir}/biom-taxonomy.tsv 

    # Finally, add the taxonomy data to your .biom file:
    biom add-metadata \
         -i {params.out_dir}/feature-table.biom \
         -o {params.out_dir}/feature-table-with-taxonomy.biom \
         --observation-metadata-fp {params.out_dir}/biom-taxonomy.tsv \
         --sc-separated taxonomy

    # Creating a TSV BIOM table
    biom convert \
           -i  {params.out_dir}/feature-table-with-taxonomy.biom  \
           -o  {params.out_dir}/feature-table-with-taxonomy.biom.tsv \
           --to-tsv
    """

shell:
    """
    set +u
    {params.conda_activate}
    set -u

    # ----- Annotate your enzymes, KOs and pathways by adding a description column ------#
    # EC
    add_descriptions.py -i {input.ec} -m EC -o {output.ec}

    # Metacyc Pathway
    add_descriptions.py -i {input.pathway} -m METACYC -o {output.pathway}

    # KO
    add_descriptions.py -i {input.ko} -m KO -o {output.ko} 

    # Unizip the metagenome contribution files - these files describe the micribes contribution the function profiles
    #find {params.outdir} -type f -name "*contrib.tsv.gz" -exec gunzip {{}} \;
    """

shell:
    """
        # Create an empty file
        mkdir -p {params.outdir} && touch {output.ko}
    """