---
# System prepended metadata

title: GreenGenes2

---

# GreenGenes2 
[Web](https://greengenes2.ucsd.edu/) info GreenGenes2
### Instalando 

Primero lo instalé en la versión más nueva de qiime. (qiime2-amplicon-2024.2) Daba un error de:

> Error: QIIME 2 has no plugin/command named 'greengenes2'.

Cree un nuevo enviromennt con una versíon más antigua de qiime 

```
wget https://data.qiime2.org/distro/core/qiime2-2023.2-py38-linux-conda.yml
conda env create -n qiime2-2023.2 --file qiime2-2023.2-py38-linux-conda.yml
```

Instalé de nuevo el plugin 

```
pip install q2-greengenes2
```

Y ahora parece que funciona :) al menos no sale el error al ejecutar 

`qiime greengenes2`  :woman-with-bunny-ears-partying: 

## Uando greengenes2 

### Datos de amplicones desde la feature-table 

Uso los datos de Nutrision la feature table 

`conda activate qiime2-2023.2 `

Descargamos la ref adecuada:

`wget http://ftp.microbio.me/greengenes_release/2022.10/2022.10.backbone.full-length.fna.qza` 


```bash
 qiime greengenes2 non-v4-16s \
    --i-table featureTable.qza \
    --i-sequences featureSeqs.qza \
    --i-backbone 2022.10.backbone.full-length.fna.qza \
    --o-mapped-table nutrision.gg2.biom.qza \
    --o-representatives nutrision.gg2.fna.qza
```

Siguiente paso: 

```bash
qiime greengenes2 taxonomy-from-table \
     --i-reference-taxonomy 2022.10.taxonomy.asv.nwk.qza \
     --i-table nutrision.gg2.biom.qza \
     --o-classification nutrision.gg2.taxonomy.qza
```

Cambiar formato para visualizar:

```
qiime metadata tabulate \
  --m-input-file nutrision.gg2.taxonomy.qza \
  --o-visualization nutrision.gg2.taxonomy.qz
```


```
qiime taxa barplot --i-table nutrision.gg2.biom.qza \
                   --i-taxonomy nutrision.gg2.taxonomy.qza \
                   --m-metadata-file metadata.tsv \
                   --o-visualization taxa_barplot.qzv
                  
``` 
                   
```
qiime feature-table summarize \
  --i-table nutrision.gg2.biom.qza \
  --o-visualization nutrision.gg2.biom.qzv \
  --m-sample-metadata-file metadata.tsv
```

Ya tengo mis datos analizados con GG2 :hamster: 

## Shotgun data 
Descargar WoRL2 

> Nota: No estoy segura si descargue la adecuada hay varias opciones para bowtie y no dice nada ._. 

Me descargué todos los archivos y los metí en una carpeta llamada `WoLr2` pero me daba un error de que no conseguía el indice. 

```bash
wget --no-check-certificate https://ftp.microbio.me/pub/wol2/databases/bowtie2/WoLr2.4.bt2l
#Aqui toda la carpeta 
wget --no-check-certificate https://ftp.microbio.me/pub/wol2/databases/bowtie2/
```

**Usar bowtie**
 
Minimize output alignment file size (no header, no unaligned, no alignment
sequences and scores, compress SAM file):

Con mis datos 

Copiar un par de datos de AI4food limpios de host

```
scp -p laura@172.16.200.189:/disk2/Grogu/victor/AI4FOOD/CLEAN_READS/V1_783_1.fastq /home/lauramarcos/Laura/gg2

scp -p laura@172.16.200.189:/disk2/Grogu/victor/AI4FOOD/CLEAN_READS/V1_783_2.fastq /home/lauramarcos/Laura/gg2
```

> Si nofunciona el indice, desacargar de `/genomes` all.fna y hacer el indice yo. 

Descargar genomas concatenados y descomprimir.

```
wget --no-check-certificate https://ftp.microbio.me/pub/wol2/genomes/all.fna.xz
mkdir -p databases/bowtie2
xz -d -k all.fna.xz
```

Construir [indice](https://www.metagenomics.wiki/tools/bowtie2/index)
```bash
bowtie2-build --seed 42 --threads 9 all.fna.xz /home/lauramarcos/Compubio/databases/bowtie2/WoLr2
#añadir al path
export BOWTIE2_INDEXES=/home/lauramarcos/Compubio/databases/bowtie2/
cp *.bt2l $BOWTIE2_INDEXES
```
* No pude por falta de espacio... 
>Sigo intentando con el pre-built sino lo hago en BIOALM y me lo paso, pero es tontería


Repito el alineamiento cn bowtie usando el indice pre-built de WoLr2 que me descargué y me di cuenta que el error puede ser que no estaba definiendo el `basename` es decir debía poner WoLr2 después del directorio, que le puse el mismo nombre :raised_hand_with_fingers_splayed: 

```bash
Este el codigo
bowtie2 -p # -x WoLr2.4.bt2l -1 input.R1.fq -2 input.R2.fq --very-sensitive \
  --no-head --no-unal | cut -f1-9 | sed 's/$/\t*\t*/' | gzip > output.sam.gz

#Mi codigo con 7 cores y todo el directorio del indice ._. 
bowtie2 -p 7 -x ~lauramarcos/Laura/gg2/WoLr2/WoLr2 -1 V1_783_1.fastq -2 V1_783_2.fastq --very-sensitive --no-head --no-unal | cut -f1-9 | sed 's/$/\t*\t*/' | gzip > output.sam.gz
```
> Si esto funciona arreglar para que este el indice exportado en el path y no tenga que poner todo el directorio y solo WoLr2/WoLr2 o algo así. Creo que así funcionaría `export BOWTIE2_INDEXES=/home/lauramarcos/Compubio/databases/bowtie2/
cp *.bt2l $BOWTIE2_INDEXES` Entonces al llamarla por WoLr2 saltaría. Probar antes de llevar al cluster. 


#### Continuación con qiime2 
Ver [aquí](https://forum.qiime2.org/t/introducing-greengenes2-2022-10/25291) info de uso con qiime2.

Instalar Woltka 

```bash
pip install woltka
#No se si hace falta 
qiime dev refresh-cache
```

Crear formato biom como imput para qiime2. Creamos **OGUs** (operational genomic unit) (Zhu et al., 2022) is a notion we proposed to define the minimum unit of microbiome composition allowed by shotgun metagenomic data. - con Woltka

Primero descomprimí el archivo autput del alinemaniento 
`gzip -d output.sam.gz`


```bash
woltka classify -i output.sam -o ogu.biom --no-demux
```
>Sino ponia --no-demux luego no me quedaba el ID de la meustra y petaba 

Poner tabla en formato qiime
```bash
conda activate qiime2-2023.2  
qiime tools import \
     --input-path ogu.biom \
     --output-path ogu.biom.qza \
     --type FeatureTable[Frequency]
```

biom convert -i ogu.biom -o feature-table.tsv --to-tsv

```
qiime metadata tabulate \
  --m-input-file ogu.biom.qza \
  --o-visualization ogu.biom.qzv
```


Filtrar features que solo estan en GG2 y luego asignar la taxonomia adecuada. 

```
qiime greengenes2 filter-features \
     --i-feature-table ogu.biom.qza \
     --i-reference 2022.10.taxonomy.asv.nwk.qza \
     --o-filtered-feature-table woltka_gg2.biom.qza

qiime greengenes2 taxonomy-from-table \
     --i-reference-taxonomy 2022.10.taxonomy.asv.nwk.qza \
     --i-table woltka_gg2.biom.qza \
     --o-classification woltka_gg2.taxonomy.qza
```

Cambiar formato para visualizar:

```
qiime metadata tabulate \
  --m-input-file woltka_gg2.taxonomy.qza \
  --o-visualization woltka_gg2.taxonomy.qzv
```

Mas visualizaciones
```
qiime taxa barplot --i-table woltka_gg2.biom.qza \
                   --i-taxonomy woltka_gg2.taxonomy.qza \
                   --m-metadata-file metadata.tsv \
                   --o-visualization taxa_barplot.qzv
                  
``` 
                  
```
qiime feature-table summarize \
  --i-table woltka_gg2.biom.qza \
  --o-visualization woltka_gg2.biom.qzv \
  --m-sample-metadata-file metadata.tsv
```


:::spoiler 
Ejemplo de uso 

redbiom fetch samples-contained --context $ctx | \
     head | \
     redbiom fetch samples \
         --context $ctx \
         --resolve-ambiguities merge \
         --output woltka.example.biom
         
biom convert -i woltka.example.biom -o woltka.example.tsv --to-tsv
       
         
qiime tools import \
     --input-path woltka.example.biom \
     --output-path woltka.example.biom.qza \
     --type FeatureTable[Frequency]
:::


---
Lanzar esto en el servidor

user: ljmarcos 

ssh ljmarcos@login2.ccc.uam.es

### 1. Pedir que instalen/descarguen

- [x] bowtie2 (ya estará)
- [x] woltka
- [x] qiime2 (ya estará), pero revisar plugin greengenes2 
- [x] Index WorL2 bowtie y pedir que lo pongan en el path
- [x] Index GG2 para qiime - creo que esto me lo paso yo misma (Pesa poco)
- [x] Metamobilepicker 


### 1.1 Pasar index bowtie 

```
scp *.bt2l ljmarcos@login2.ccc.uam.es:/home/ljmarcos/WoLr2
```
### 2. Pasar datos filtrados de host de AI4Food - sólo V1 

```
scp *.fastq ljmarcos@login2.ccc.uam.es:/home/ljmarcos/clean_data
```

### 3. Preprar scripts 

Deberia hacer uno para SLURM y luego uno para bowtie, otro para woltka y otro para qiime2. 

#### Bowtie2

`bowtie2 -p 7 -x ~lauramarcos/Laura/gg2/WoLr2/WoLr2 -1 *1.fastq -2 *_2.fastq --very-sensitive --no-head --no-unal | cut -f1-9 | sed 's/$/\t*\t*/' | gzip > output.sam.gz`

#### Woltka 

```
woltka classify -i *.sam -o ogu.biom --no-demux
```

#### Qiime2 
Sale mas a cuenta hacer esto en Grogu, es rapido. 

:::spoiler
Sería así
```bash=1
#!/bin/bash
#SBATCH -p bioinfo                      #partition/queue name
#SBATCH --job-name=woltka                #Job name
#SBATCH -N 1                   #Un nodo
#SBATCH -n 40                   #Cuarenta cores
#SBATCH -t 05:00:00             #Time limit hrs:min:sec
#SBATCH --output=/home/ljmarcos/qiime/log-%j.o #Log de salida
#SBATCH --error=/home/ljmarcos/qiime/log-%j.e #Log de errores
#SBATCH --mail-user=judith.marcos@alimentacion.imdea.org
#SBATCH --mail-type=END,FAIL

module load qiime/2023.2
conda activate qiime2-2023.2  
qiime tools import \
     --input-path ogu.biom \
     --output-path ogu.biom.qza \
     --type FeatureTable[Frequency]

#Filtrar features que solo estan en GG2 y luego asignar la taxonomia adecuada. 

qiime greengenes2 filter-features \
     --i-feature-table ogu.biom.qza \
     --i-reference 2022.10.taxonomy.asv.nwk.qza \
     --o-filtered-feature-table woltka_gg2.biom.qza

qiime greengenes2 taxonomy-from-table \
     --i-reference-taxonomy 2022.10.taxonomy.asv.nwk.qza \
     --i-table woltka_gg2.biom.qza \
     --o-classification woltka_gg2.taxonomy.qza

#Cambiar formato para visualizar:

qiime metadata tabulate \
  --m-input-file woltka_gg2.taxonomy.qza \
  --o-visualization woltka_gg2.taxonomy.qzv

```
:::

### 4. Lanzar trabajos 

Script bowtie2:

```perl=1
#!/usr/bin/perl
foreach $ar (@ARGV)
{
    $r1 = $ar;
    $r2 = $ar;
    $r2 =~ s/_1/_2/;

    $name = $ar;
    $name =~ s/_1.fastq//;

    system("bowtie2 -p 20 -x /usr/local/BBDD/WoLr2/WoLr2 -1 $r1 -2 $r2 --very-sensitive --no-head --no-unal | cut -f1-9 | sed 's/$/\t*\t*/' | gzip > $name.sam.gz");

}
```

Script Slurm 

```bash=1
#!/bin/bash
#SBATCH -p bioinfo                      #partition/queue name
#SBATCH --job-name=alineamiento_WoRL2                #Job name
#SBATCH -N 1                   #Un nodo
#SBATCH -n 20                   #Veinte cores
#SBATCH -t 2000:00:00             #Time limit hrs:min:sec
#SBATCH --initialdir="/home/proyectos/imdeaalim/ljmarcos/clean_data"
#SBATCH --output=/home/ljmarcos/WoRL2/log-%j.o #Log de salida
#SBATCH --error=/home/ljmarcos/WoRL2/log-%j.e #Log de errores
#SBATCH --mail-user=judith.marcos@imdea.org
#SBATCH --mail-type=END,FAIL

module load bowtie2/2.5.3
perl bowtie2.pl *_1.fastq

echo "Job" $SLURM_JOB_NAME 
echo "ID" $SLURM_JOBID
exit 0
```

**Ejecutar** :computer: 

Conectar al servidor y ejecutar script slurm :snail:

```bash
ssh ljmarcos@login2.ccc.uam.es
cd .. 
cd proyectos/imdeaalim/ljmarcos/clean_data/
sbatch -A imdeaalim_serv -p bioinfo slurm_bowtie2.sh
```

## Woltka

Para sacar la lista de nombres de muestras: 

```bash=1
#!/bin/bash
# Directorio de entrada
directorio="."

# Archivo de salida
archivo_salida="samples.txt"

# Iterar sobre los archivos en el directorio
for archivo in "$directorio"/*; do
    # Extraer solo el nombre del archivo sin la extensión
    nombre_archivo=$(basename "$archivo" | cut -d. -f1)
    # Concatenar los nombres con una coma
    nombres="$nombres$nombre_archivo,"
done

# Eliminar la última coma
nombres=${nombres%,}

# Escribir los nombres de archivo en el archivo de salida
echo "$nombres" > "$archivo_salida"

echo "Proceso completado. Se ha generado el archivo '$archivo_salida'."
```

Script de slurm :snail: 

```bash=1
#!/bin/bash
#SBATCH -p bioinfo                      #partition/queue name
#SBATCH --job-name=alineamiento_WoRL2                #Job name
#SBATCH -N 1                   #Un nodo
#SBATCH -n 20                   #Veinte cores
#SBATCH -t 2000:00:00             #Time limit hrs:min:sec
#SBATCH --initialdir="/home/proyectos/imdeaalim/ljmarcos/clean_data"
#SBATCH --output=/home/ljmarcos/woltka/log-%j.o #Log de salida
#SBATCH --error=/home/ljmarcos/woltka/log-%j.e #Log de errores
#SBATCH --mail-user=judith.marcos@imdea.org
#SBATCH --mail-type=END,FAIL

module load woltka/0.1.6
woltka classify -i /home/proyectos/imdeaalim/ljmarcos/clean_data/alignment -s samples.txt -o ogu.biom --no-demux
```

Conectar al servidor y ejecutar script slurm :snail:

```bash
ssh ljmarcos@login2.ccc.uam.es
cd .. 
cd proyectos/imdeaalim/ljmarcos/clean_data/
sbatch -A imdeaalim_serv -p bioinfo slurm_woltka.sh
```
Me da un error y aunque diga que admite archivos comprimidos, ya me pasó en la prueba que no era así... modifico el script para que primero descomprima todo. 

```bash=1
#!/bin/bash
#SBATCH -p bioinfo                      #partition/queue name
#SBATCH --job-name=woltka                #Job name
#SBATCH -N 1                   #Un nodo
#SBATCH -n 40                   #Veinte cores
#SBATCH -t 72:00:00             #Time limit hrs:min:sec
#SBATCH --output=/home/ljmarcos/woltka/log-%j.o #Log de salida
#SBATCH --error=/home/ljmarcos/woltka/log-%j.e #Log de errores
#SBATCH --mail-user=judith.marcos@imdea.org
#SBATCH --mail-type=END,FAIL

#Cargar woltka
module load woltka/0.1.6

# Directorio donde se encuentran los archivos de alineamiento
directorio_alignment="/home/proyectos/imdeaalim/ljmarcos/clean_data/alignment"

# Descomprimir archivos en el directorio de alineamiento
gunzip -r "$directorio_alignment"

# Ejecutar el comando Woltka
woltka classify -i "$directorio_alignment" -s samples.txt -o ogu.biom --no-demux
```
Me sale un error raro, repito así: 



```bash=1
#!/bin/bash
#SBATCH -p bioinfo                      #partition/queue name
#SBATCH --job-name=woltka                #Job name
#SBATCH -N 1                   #Un nodo
#SBATCH -n 40                   #Cuarenta cores
#SBATCH -t 72:00:00             #Time limit hrs:min:sec
#SBATCH --output=/home/ljmarcos/woltka/log-%j.o #Log de salida
#SBATCH --error=/home/ljmarcos/woltka/log-%j.e #Log de errores
#SBATCH --mail-user=judith.marcos@imdea.org
#SBATCH --mail-type=END,FAIL

#Cargar woltka
module load woltka/0.1.6

# Ejecutar el comando Woltka
woltka classify -i /alignment -s samples.txt -o ogu.biom --no-demux
```


:memo: Actulización 20 de mayo 2024 
El error era porque algunos archivos no se habian alineado, después demucho tiempo perdido, me dicuenta que habíaun error en el script de perl al poner que sustituyera *_1* por *_2* me cambiaba el nomnre de l amuestra y no existia y hacia el alinemaiento pero con el archivo vacio... Me dic centa que todas  la smuestras llamadas "V1_1XX" estaban mal (Total 29 muestras) y esto l oque hacía es que lueg woltka no corriera. Estoy repitiendo con el script corregido (En vez de *_1* pone *_1.fastq") a ver si así funciiiona. 
La prueba con 4 alineamientos buenos funcionó, así que tengo fe. :pray: 

---

Me escriben del cluster pidiendo que modifique el script tal que así: 

```
# Cargar woltka.
module load woltka/0.1.6

# Copiar a temporal los ficheros de entrada necesarios para woltka.
mkdir  /temporal/$SLURM_JOB_USER/$SLURM_JOB_ID
cp -r alineamientos  /temporal/$SLURM_JOB_USER/$SLURM_JOB_ID
cd   /temporal/$SLURM_JOB_USER/$SLURM_JOB_ID

# Ejecutar el comando Woltka.
woltka classify -i alineamientos -f sam -o alineamientos/resultados/ogu.biom

# Mover a imdeaalim los resultados del cálculo.
mv alineamientos/resultados/ogu.biom /lustre/NodoBIO/imdeaalim/ljmarcos/clean_data/alineamientos/resultados

# Eliminar los archivos temporales.
rm -rf /temporal/$SLURM_JOB_USER/$SLURM_JOB_ID
```

El correo en cuestión: 
:::spoiler

Estimada Laura Judith:
    Hemos comprobado que sus trabajos realizan un alto impacto en el sistema de ficheros centralizado, por lo que le rogamos realice una ligera modificación en sus scripts. De hecho es muy posible que con esta modificación encuentre una mejora de rendimiento en sus trabajos.
    Centrándonos en el trabajo que hemos tenido que parar por el comportamiento que hemos indicado, en su script tiene tras las líneas para SLURM el siguiente código:
```
#Cargar woltka
module load woltka/0.1.6

# Ejecutar el comando Woltka
woltka classify -i alineamientos -f sam -o alineamientos/resultados/ogu.biom
```
Necesitaríamos que lo cambiara por:
```
# Cargar woltka.
module load woltka/0.1.6

# Copiar a temporal los ficheros de entrada necesarios para woltka.
mkdir  /temporal/$SLURM_JOB_USER/$SLURM_JOB_ID
cp -r alineamientos  /temporal/$SLURM_JOB_USER/$SLURM_JOB_ID
cd   /temporal/$SLURM_JOB_USER/$SLURM_JOB_ID

# Ejecutar el comando Woltka.
woltka classify -i alineamientos -f sam -o alineamientos/resultados/ogu.biom

# Mover a imdeaalim los resultados del cálculo.
mv alineamientos/resultados/ogu.biom /lustre/NodoBIO/imdeaalim/ljmarcos/clean_data/alineamientos/resultados

# Eliminar los archivos temporales.
rm -rf /temporal/$SLURM_JOB_USER/$SLURM_JOB_ID
```
Si se fija lo que intentamos hacer es que calcule en /temporal, que es un directorio local de las máquinas de cálculo accesible únicamente en el momento del cálculo, pero al ser local a la máquina ofrece mayor rendimiento.
   Para ello creamos un directorio propio en este sistema, si bien Slurm lo suele crear por defecto, y copiamos lo que necesitamos en él. Una vez finalizada la copia nos movemos a ese directorio y se ejecuta woltka.
   Al final se mueven los resultados al directorio original.
   Tenga presente que apenas conocemos el funcionamiento de woltka, por lo que por favor modifique estas líneas para copiar a temporal otros ficheros si son necesarios y para luego mover a su directorio de NodoBIO lo que considere que necesita.
    Por supuesto si necesita cualquier ayuda en este sentido, indíquenoslo.
    Lamentamos las molestias ocasionadas.
    Reciba un cordial saludo,
                                Pablo Sanz

:::

# Actualización Qiime/GreenGenes2 

Sería esta BD: `2024.09.backbone.full-length.fna.qza ` Esta [aquí](http://ftp.microbio.me/greengenes_release/current/).
Y la versión de qiime2 para amplicones `qiime2-amplicon-2024.2` 

```
# Ejecutar qiime2 

conda activate qiime2-amplicon-2024.2 

qiime tools import \
  --type 'SampleData[PairedEndSequencesWithQuality]' \
  --input-path manifiest.csv \
  --output-path demux.qza \
  --input-format PairedEndFastqManifestPhred33 \

  qiime demux summarize \
  --i-data demux.qza \
  --o-visualization demux.qzv

# DADA2
#screen qiime
  qiime dada2 denoise-paired \
  --i-demultiplexed-seqs demux.qza \
  --p-n-threads 10 \
  --p-trim-left-f 25 \
  --p-trim-left-r 20 \
  --p-trunc-len-f 250 \
  --p-trunc-len-r 235 \
  --o-table table.qza \
  --o-representative-sequences rep-seqs.qza \
  --o-denoising-stats denoising-stats.qza \
  --p-min-fold-parent-over-abundance 10 
  

#Con ptrunc-len-r a 240 mejora elporcentaje de merged de 30 a 50% todavia hay muchas quimeras, 
pruebo con 210 y luego veré con el parametro de filtrar quimeras, tambien quité mas del trim-f dejandolo en 25 
210 mal, a ver 235 sino 240 y subimos el parent-over a 8 

#Visualizaciones
qiime feature-table summarize \
  --i-table table.qza \
  --o-visualization table.qzv \
  --m-sample-metadata-file metadata.tsv

  qiime feature-table tabulate-seqs \
  --i-data rep-seqs.qza \
  --o-visualization rep-seqs.qzv

qiime metadata tabulate \
  --m-input-file denoising-stats.qza \
  --o-visualization denoising-stats.qzv


# GreenGenes2 para taxonomía 
#Instalé el plugin de GreenGenes2 para qiime 

pip install q2-greengenes2

#Usé el nuevo backbone actualizado además es para qiime-amplicon-2024 :) 

qiime greengenes2 non-v4-16s \
    --i-table table.qza \
    --i-sequences rep-seqs.qza \
    --i-backbone 2024.09.backbone.full-length.fna.qza \
    --o-mapped-table results.gg2.biom.qza \
    --o-representatives results.gg2.fna.qza

qiime greengenes2 taxonomy-from-table \
     --i-reference-taxonomy 2024.09.taxonomy.asv.nwk.qza \
     --i-table results.gg2.biom.qza \
     --o-classification results.gg2.taxonomy.qza

     qiime metadata tabulate \
  --m-input-file results.gg2.taxonomy.qza \
  --o-visualization results.gg2.taxonomy.qzv

qiime taxa barplot --i-table results.gg2.biom.qza \
                   --i-taxonomy results.gg2.taxonomy.qza \
                   --m-metadata-file metadata.tsv \
                   --o-visualization taxa_barplot.qzv

qiime feature-table summarize \
  --i-table results.gg2.biom.qza \
  --o-visualization results.gg2.biom.qzv \
  --m-sample-metadata-file metadata.tsv

  # Colapsar a nivel de género 

#Comandos en qiime para colapsar las features 16S 
qiime taxa collapse \
 --i-table results.gg2.biom.qza \
 --i-taxonomy results.gg2.taxonomy.qza \
 --p-level 6 \
 --o-collapsed-table 16Scollapsed.qza
```