> For the complete documentation index, see [llms.txt](https://jodyphelan.gitbook.io/tb-profiler-pt/llms.txt). Markdown versions of documentation pages are available by appending `.md` to page URLs; this page is available as [Markdown](https://jodyphelan.gitbook.io/tb-profiler-pt/development/escrevendo-um-script-collate-personalizado.md).

# Escrevendo um script collate personalizado

Os arquivos resultantes de cada execução individual são armazenados por padrão em arquivos json (terminando em **.results.json**). Estes podem ser agrupados em um arquivo mesclado para exibir as mutações de resistência às drogas usando o `tb-profiler collate`. No entanto, você pode querer personalizar quais informações são extraídas dos arquivos resultantes e exibidas no output agrupado (collated). Esta página fornecerá um esboço de como escrever seu próprio analisador.

## Esqueleto do script

Abra seu editor favorito, cole o código a seguir e salve o arquivo como **tbprofiler\_custom\_script.py**

```python
#! /usr/bin/env python

import sys
import csv
import json
import argparse
import os
from collections import defaultdict
from tqdm import tqdm
import pathogenprofiler as pp
import tbprofiler as tbprofiler

def main(args):
    bed_file = "%s/share/tbprofiler/%s.bed" % (sys.base_prefix,args.db)
    locus_tag2drugs = tbprofiler.get_lt2drugs(bed_file)

    if args.samples:
        samples = [x.rstrip() for x in open(args.samples).readlines()]
    else:
        samples = [x.replace(args.suffix,"") for x in os.listdir(args.dir) if x[-len(args.suffix):]==args.suffix]

    for s in tqdm(samples):
        data = json.load(open(pp.filecheck("%s/%s%s" % (args.dir,s,args.suffix))))


parser = argparse.ArgumentParser(description='tbprofiler script',formatter_class=argparse.ArgumentDefaultsHelpFormatter)
parser.add_argument('--samples',type=str,help='File with samples')
parser.add_argument('--dir',default="results/",type=str,help='Directory containing results')
parser.add_argument('--db',default="tbdb",type=str,help='Database name')
parser.add_argument('--suffix',default=".results.json",type=str,help='File suffix')
parser.set_defaults(func=main)

args = parser.parse_args()
args.func(args)
```

Isso formará a estrutura básica do analisador. Vamos ver o que cada seção faz.

**Linhas 3-11**: Aqui importamos bibliotecas que são úteis no script.

**Linhas 13**: esta linha inicia a definição da função principal (`main` function). Recebe um argumento na forma de um objeto de argumentos que é construído usando o analisador de argumentos.

**Linha 14**: A localização do arquivo bed  é armazenada na variável `bed_file`. Isso é construído usando o prefixo do banco de dados fornecido nos argumentos da linha de comando (tbdb por padrão).

**Linha 15**: Usamos a função `tbprofiler.get_lt2drugs` para extrair o locus\_tag para o mapeamento de droga associado em um dicionário com o seguinte formato:

```python
{
    "Rv0667":["rifampicin"],
    "Rv1484":["isoniazid","ethionamide"],
}
```

**Linha 17-20**: Se o argumento `--samples` foi fornecido, o arquivo é lido e armazenado como uma lista. Caso contrário, a lista é derivada da procura do diretório de resultados (especificado com `--dir`).

**Linhas 22-23**: Um loop é construído para percorrer as amostras. Para cada iteração, o arquivo resultante json é carregado em uma estrutura de dicionário. A estrutura de valor-chave é explicada abaixo.

**Linhas 26-34**: Essas linhas analisam os argumentos fornecidos ao script (por exemplo -`--samples`). Você pode adicionar novos argumentos e eles aparecerão como propriedades do objeto `args` na função principal (`main` function).

#### Estrutura de dados do resultado

Os arquivos de resultados individuais são lidos em um dicionário na linha 23. Aqui estão algumas (mas não todas) as chaves interessantes no dicionário:

```python
{
    "main_lin": "lineage4",
    "sublin": "lineage4.9",
    "dr_variants": [
      {
        "sample": "por1A",
        "gene_name": "rpoB",
        "chr": "Chromosome",
        "genome_pos": 761155,
        "type": "missense",
        "change": "p.Ser450Leu",
        "freq": 1,
        "nucleotide_change": "761155C>T",
        "locus_tag": "Rv0667",
        "gene": "rpoB",
        "_internal_change": "450S>450L",
        "drug": "rifampicin",
        "literature": "10.1128/AAC.01093-18",
        "confidence": "high"
      }
    ],
    "other_variants": [
      {
        "sample": "DRR034528",
        "gene_name": "gyrA",
        "chr": "Chromosome",
        "genome_pos": 7362,
        "type": "missense",
        "change": "p.Glu21Gln",
        "freq": 1,
        "nucleotide_change": "7362G>C",
        "locus_tag": "Rv0006",
        "gene": "gyrA",
        "_internal_change": "21E>21Q"
      }
    ]
}
```

As mutações são divididas em duas listas “**dr\_variants”** e “**other\_variants**”. As “**dr\_variants”** estão todas associadas à resistência às drogas. As outras variantes (“**other\_variants”**) são mutações em genes resistentes às drogas, mas que não foram (ainda!) associadas à resistência. Algumas das outras mutações podem ser mutações filogenéticas que não têm nada a ver com resistência, enquanto outras podem ser, mas podem não ter evidências suficientes ainda.

No momento, se você executar o script e apontar para um diretório com alguns arquivos de resultado usando `--dir` ele lerá todos os dados, mas não fará nada com eles.

## **Roteiro de exemplo de trabalho**

Como exemplo, estenderei sua funcionalidade para gerar um arquivo CSV simples com o nome da amostra e o tipo de resistência à droga.

Primeiro, criarei um objeto de gravação CSV. Adicione o seguinte código antes do for loop.

```python
OUT = open(args.out,"w")
writer = csv.DictWriter(OUT, fieldnames = ["sample","dr-class"])
writer.writeheader()
```

Feche o identificador de arquivo após o término do loop usando

```python
OUT.close()
```

Usamos a propriedade `args.out` ao criar o filehandle. Para que isso funcione, temos que adicionar um argumento extra ao analisador de argumentos. Adicione o seguinte para o código junto com os outros argumentos:

```python
parser.add_argument('--out', type=str, help='Name of CSV output', required=True)
```

Agora vamos adicionar naquele código que determina a classe de resistência às drogas. Por padrão, o tb-profiler classifica apenas as amostras em: sensível, MDR, XDR e droga resistente. Vamos estender isso para encontrar pré-MDR e pré-XDR.

Primeiro, podemos criar um [conjunto](https://docs.python.org/3.7/library/stdtypes.html#set-types-set-frozenset) contendo uma lista dos medicamentos aos quais a amostra é resistente. Podemos fazer isso usando uma [compreensão de lista](https://docs.python.org/3.7/howto/functional.html?highlight=list%20comprehension#generator-expressions-and-list-comprehensions) que é convertida em um conjunto. Adicione o seguinte código dentro do for loop:

```python
resistant_drugs = set([d["drug"] for d in data["dr_variants"]])
```

#### Definições das classes de drogas

As amostras podem ser classificadas em diferentes tipos usando as seguintes definições

| Tipo     | Droga resistência                                                            |
| -------- | ---------------------------------------------------------------------------- |
| Sensível | Sem droga resistance                                                         |
| Pré-MDR  | Rifampicina **ou** isoniazida                                                |
| MDR      | Rifampicina **e** isoniazida                                                 |
| Pré-XDR  | MDR **e** (qualquer fluoroquinolona **ou** qualquer segunda-linha injetável) |
| XDR      | MDR **e** (qualquer fluoroquinolona e qualquer segunda-linha injetável)      |
| Outras   | Resistência para qualquer droga mas nenhuma das categorias acima             |

Codificando isso

Podemos ver imediatamente que as definições são apenas algumas operações booleanas. Podemos armazenar os valores booleanos procurando as drogas no conjunto `resistant_drugs` que criamos antes. Primeiro, criaremos uma lista de fluoroquinolonas (FLQ) e de segunda linha injetáveis (SLI). Adicione este código em qualquer lugar na função principal, mas antes do loop:

```python
FLQ_set = set(["moxifloxacin","levofloxacin","ciprofloxacin","ofloxacin"])
SLI_set = set(["amikacin","capreomycin","kanamycin"])
```

Em seguida, crie os valores booleanos. Coloque-os dentro do loop depois de criar a variável de `data`:

```python
resistant_drugs = set([d["drug"] for d in data["dr_variants"]])
rif = "rifampicin" in resistant_drugs
inh = "isoniazid" in resistant_drugs
flq = len(FLQ_set.intersection(resistant_drugs)) > 0
sli = len(SLI_set.intersection(resistant_drugs)) > 0
```

Para testar a rifampicina e a isoniazida, verificamos simplesmente se elas existem no `resistant_drugs`. Para verificar a resistência FLQ e SLI, encontramos as interseções entre `FLQ_set` e `resistant_drugs` e count. Se o valor for maior que 0, então temos resistência nessa categoria.

Agora estamos prontos para realizar a construção das classes usando nossas variáveis booleanas. Podemos percorrer a tabela acima usando uma série de instruções if-else.

```python
if len(resistant_drugs)==0:
    drtype = "Sensitive"
elif (rif and not inh) or (inh and not rif):
    drtype = "Pre-MDR"
elif (rif and inh) and (not flq and not sli):
    drtype = "MDR"
elif (rif and inh) and ( (flq and not sli) or (sli and not flq) ):
    drtype = "Pre-XDR"
elif (rif and inh) and (flq and sli):
    drtype = "XDR"
else:
    drtype = "Other"
```

Armazenamos a classe de resistência aos medicamentos na variável `drtype`. Agora, tudo o que resta fazer é gravar isso no arquivo CSV. Podemos fazer isso usando:

```python
writer.writerow({"sample":s, "dr-class":drtype})
```

Seu script agora deve ter a seguinte aparência:

```python
 #! /usr/bin/env python

import sys
import csv
import json
import argparse
import os
from collections import defaultdict
from tqdm import tqdm
import pathogenprofiler as pp
import tbprofiler as tbprofiler

def main(args):
    bed_file = "%s/share/tbprofiler/%s.bed" % (sys.base_prefix,args.db)
    locus_tag2drugs = tbprofiler.get_lt2drugs(bed_file)

    if args.samples:
        samples = [x.rstrip() for x in open(args.samples).readlines()]
    else:
        samples = [x.replace(args.suffix,"") for x in os.listdir(args.dir) if x[-len(args.suffix):]==args.suffix]

    FLQ_set = set(["moxifloxacin","levofloxacin","ciprofloxacin","ofloxacin"])
    SLI_set = set(["amikacin","capreomycin","kanamycin"])

    OUT = open(args.out,"w")
    writer = csv.DictWriter(OUT, fieldnames = ["sample","dr-class"])
    writer.writeheader()

    for s in tqdm(samples):
        data = json.load(open(pp.filecheck("%s/%s%s" % (args.dir,s,args.suffix))))
        resistant_drugs = set([d["drug"] for d in data["dr_variants"]])

        resistant_drugs = set([d["drug"] for d in data["dr_variants"]])
        rif = "rifampicin" in resistant_drugs
        inh = "isoniazid" in resistant_drugs
        flq = len(FLQ_set.intersection(resistant_drugs)) > 0
        sli = len(SLI_set.intersection(resistant_drugs)) > 0

        if len(resistant_drugs)==0:
            drtype = "Sensitive"
        elif (rif and not inh) or (inh and not rif):
            drtype = "Pre-MDR"
        elif (rif and inh) and (not flq and not sli):
            drtype = "MDR"
        elif (rif and inh) and ( (flq and not sli) or (sli and not flq) ):
            drtype = "Pre-XDR"
        elif (rif and inh) and (flq and sli):
            drtype = "XDR"
        else:
            drtype = "Other"

        writer.writerow({"sample":s, "dr-class":drtype})

    OUT.close()

parser = argparse.ArgumentParser(description='tbprofiler script',formatter_class=argparse.ArgumentDefaultsHelpFormatter)
parser.add_argument('--out', type=str, help='Name of CSV output', required=True)
parser.add_argument('--samples', type=str, help='File with samples')
parser.add_argument('--dir', default="results/", type=str, help='Directory containing results')
parser.add_argument('--db', default="tbdb", type=str, help='Database name')
parser.add_argument('--suffix', default=".results.json", type=str, help='File suffix')
parser.set_defaults(func=main)

args = parser.parse_args()
args.func(args)
```

Tente executá-lo e verifique se o output parece como deveria:

```bash
python tbprofiler_custom_script.py --dir /path/to/dir --out test.csv
```

Se estiver funcionando conforme o esperado, tente adicionar outras informações ao output. Por exemplo, você pode extrair a sublinhagem acessando `data["sublin"]`. Dê uma olhada na estrutura do arquivo de resultado (você pode abri-lo com o [Firefox](https://www.mozilla.org/en-GB/firefox/new/) para direcioná-lo a um formato mais legível). Esperançosamente, você pode usar este modelo como um padrão para seus próprios usos. Existem alguns exemplos disponíveis [aqui](https://github.com/jodyphelan/tbprofiler_scripts) se você estiver interessado. Se você tiver alguma dúvida ou se tiver um script legal que gostaria de disponibilizar para a comunidade, não hesite em me contatar por e-mail ou pelo [GitHub](https://github.com/jodyphelan/TBProfiler/issues). Obrigado por ler!
