Aula – Montagem de genomas

 

1. Abra um terminal em sua máquina Linux e crie uma pasta mkdir aula_montagem e entre nela cd aula_montagem

 

2. Agora abra uma conexão com a pinguim ssh bioufmg@143.107.223.182

 

3. Entre em sua pasta cd eusoujacu crie um diretório mkdir aula_montagem e entre nela cd aula_montagem

 

4. Copie o material para essa pasta cp /home/treinamento/velvet_aula/* . (olha o ponto!)

 

5. Dá um ls e veja que você tem três arquivos, dois arquivos de sequenciamento (reads_forward.fastq e reads_reverse.fastq) e um arquivo com o genoma de referência (genoma.fasta)

 

6. Primeiro vamos olhar a qualidade dos reads com o programa fastQC, ele vai gerar várias saídas. Para ver as saídas temos uma novidade. Na biodados (http://biodados.icb.ufmg.br link you@pinguim) há um link para public_html, que mostra na web o conteúdo de todas as pastas em bioufmg. Assim, toda figura gerada lá pode ser vista com o firefox!

 

7. Para observar a qualidade das bases nos arquivos reads_forward.fastq e reads_reverse.fastq, crie o diretório sem_trim e execute o script fastqc (com –t especifique rodar com um único core/thread)

mkdir sem_trim  (o programa exige que o diretório de output seja criado antes...)

fastqc -o sem_trim -t 1 reads_forward.fastq  reads_reverse.fastq

 

8. Veja a saída abrindo a public_html, abrindo sua pasta e seguindo o caminho: aula_montagem/sem_trim/reads_forward_fastqc/fastqc_report.html para informações do reads forward e aula_montagem/sem_trim/reads_reverse_fastqc/fastqc_report.html para informações do reads reverse

 

9. Vamos filtrar esses reads com Trimmomatic assim:

/usr/local/bin/trimmomatic PE -phred33 reads_forward.fastq reads_reverse.fastq trimmed_forward_paired.fastq trimmed_forward_unpaired.fastq trimmed_reverse_paired.fastq trimmed_reverse_unpaired.fastq LEADING:3 TRAILING:3 SLIDINGWINDOW:4:15 MINLEN:36

Os parâmetros para trimming são:

LEADING: Remove bases com baixa qualidade no início da sequencia (qualidade abaixo de 3)

TRAILING: Remove bases com baixa qualidade no final da sequencia (qualidade abaixo de 3)

SLIDINGWINDOW: Utiliza uma janela deslizando com 4 bases e corta quando a qualidade média é menor que 15

MINLEN: Remove sequencias menores que 36 bases

 

10. Agora refaça o processo do fastQC com os reads filtrados

mkdir com_trim  

fastqc -o sem_trim -t 1 trimmed_forward_paired.fastq trimmed_reverse_paired.fastq

Veja o report do programa no firefox, em public_html abra sua pasta e siga o caminho:

aula_montagem/sem_trim/trimmed_forward_paired_fastqc/fastqc_report.html

aula_montagem/sem_trim/trimmed_reverse_paired_fastqc/fastqc_report.html

Melhorou?

 

11. Com os reads trimados fazer a montagem sem referencia (do zero) com o programa Velvet

 

12. Rode o velvet assim:

Crie os índices (o parâmetro numérico refere-se ao tamanho do k-mer utilizado, sempre números ímpares)

velveth montagem 31 -fastq -shortPaired -separate trimmed_forward_paired.fastq trimmed_reverse_paired.fastq

Execute a montagem

velvetg montagem/ -exp_cov auto -scaffolding yes

 

Após a execução, o Velvet apresenta algumas informações da montagem. Os nodes são o número de contigs gerados. O n50 diz que metade do genoma está representada em contigs maiores que esse valor.

 

13. Verifique o arquivo de montagem

cd montagem

more contigs.fa

 

 

14. Agora vamos fazer o outro tipo de montagem que é "com genoma de referencia" usando o programa Bowtie. Primeiro volte ao diretório aula_montagem cd ..

Rode o bowtie assim:

Crie um bowtie index usando o genoma de ancoragem (cts.fasta) fornecido:

bowtie2-build genoma.fasta genoma.build

Alinhe os reads no index do genoma de referência (genoma.build), criando um arquivo de alinhamentos reads_mapeados.sam no formato chamado SAM (-S manda criar o SAM)

bowtie2 -x genoma.build -1 trimmed_forward_paired.fastq -2 trimmed_reverse_paired.fastq -S reads_mapeados.sam

 

15. Para ver o resultado teremos que abrir um programa Java na sua máquina Linux (ou na windows) chamado Tablet. Você pode rodar ele da web usando este link. Ele vai abrir um ambiente gráfico. Se não funcionar na sua máquina pessoal, você pode fazer o download do programa aqui.

 

16. Temos que pegar o arquivo de saída do Bowtie (reads_mapeados.sam) e o genoma de referência (genoma.fasta) e trazer para a sua máquina. No Linux da aula, o melhor a fazer é achar o resultado usando o you@pinguim (link na biodados) e "copiar link", depois, no terminal da sua máquina, dar um wget, assim:

wget http://143.107.223.182/public_html/eusoujacu/aula_montagem/genoma.fasta

wget http://143.107.223.182/public_html/eusoujacu/aula_montagem/reads_mapeados.sam

 

17. Dá um ls pra ver o arquivo. E quem está no Windows? nesse caso você já sabe, é só baixar tudo como de costume e abrir pelo modo gráfico do Tablet.

 

18. Aberto o genoma de referencia e o arquivo de saída do Bowtie, você vai ver um monte de reads ancorados na referência.

Clique em: "Open Assembly" 

Na primeira janela selecione o arquivo: “reads_mapeados.sam”   

Na segunda janela selecione o genoma de referência: “genoma.fasta” 

Clique em open.

 

19. Não é aula de transcriptômica, então vamos aprender a exportar os resultados do Bowtie usando o SAMTools. De volta na pinguim:

Converta o arquivo SAM para a versão binária BAM

samtools view -S reads_mapeados.sam -b -o reads_mapeados.bam

Ordene os reads

samtools sort reads_mapeados.bam reads_mapeados_ordenados

Crie um índice para o BAM

samtools index reads_mapeados_ordenados.bam

Crie os contigs com o consenso do mapeamento

samtools mpileup -uf genoma.fasta reads_mapeados_ordenados.bam | bcftools view -cg - | vcfutils.pl vcf2fq > montagem_por_referencia.fastq

Visualize a montagem

more montagem_por_referencia.fastq