##Criar banco de dados para o bowtie

bowtie-build genoma.fasta nome_da_database

##Mapear

bowtie -n 3 -p 3 -m 1 --best --sam -1 read_1.fastq -2 read_2.fastq nome_da_database > mapped.sam

-n -> numero de mismatches aceitos 

-p -> processadores 

-m -> numero de alinhamentos por read 

--best -> mostra só o melhor alinhamento 

--sam -> output em formato SAM

##Criar formato BAM

samtools faidx genoma.fasta

samtools import genoma.fasta.fai mapped.sam mapped.bam

samtools sort mapped.bam mapped.sorted

samtools index mapped.sorted.bam 

##Quantificando

bedtools multicov -bams mapped1.sorted.bam mapped2.sorted.bam mapped3.sorted.bam ... mappedN.sorted.bam -bed genoma.gff > counts

cat aula.count | cut -f9- > data.table

##Estatística

R

##Carregar a bibliotecas

library(DESeq2)

library(gplots)

##Ler a tabela de contagens

data = read.table(file=‘data.table’, head=F, sep="\t", row.names=1)

##Definir o nome das colunas (nome das amostras)

colnames(data) = c("CI", "CII", "CIII", "CMI", "CMII", "CMIII", "RI", "RII", "RIII", "RMI", "RMIII")

##Montar uma tabela com os metadados

treatment = c("Control", "Control", "Control", "Control", "Control", "Control", "Resist", "Resist", "Resist", "Resist", "Resist")

antibiotc = c("no", "no", "no", "yes", "yes", "yes", "no", "no", "no", "yes", "yes")

metadata = data.frame(treatment=as.factor(treatment), antibiotc=as.factor(antibiotc))

rownames(metadata)=c("CI", "CII", "CIII", "CMI", "CMII", "CMIII", "RI", "RII", "RIII", "RMI", "RMIII")

metadata$treatment=relevel(metadata$treatment, "Control")

metadata$antibiotc=relevel(metadata$antibiotc, “no”)

##Fazer a análise estatística

ddsData = DESeqDataSetFromMatrix(data, colData, formula (~ treatment))

dds = DESeq(ddsData)

res = results(dds)

select = which(res$padj<0.01)

rld <- rlogTransformation(dds, blind=TRUE)

vsd <- varianceStabilizingTransformation(dds, blind=TRUE)

plotMA(res, main="DESeq2", ylim=c(-3,3), alpha=0.01)

abline(h=-1, col="blue")

abline(h=1, col="blue")

distsRL <- dist(t(assay(rld)))

mat <- as.matrix(distsRL)

rownames(mat) <- colnames(mat) <- with(colData(dds), colnames(data))

heatmap.2(mat, trace="none", margin=c(5, 5))