Delete mapeo.R

-# ENFERMEDADES MONOGENICAS
-# EQUIPO 5
-
-# Este script tiene la finalidad de obtener la anotacion en el genoma de Homo sapiens
-# de aquellas coordenadas que corresponden a los genes de interes.
-
-
-#se mandan a llamar a todas las librerias a utilizar 
-library(Biobase)
-library(IRanges)
-library(rtracklayer)
-library(GenomicRanges)
-library(Rsamtools)
-
-
-#Obtenemos la anotación de un genoma
-Para tener con que comparar el genoma que vamos a analizar, es necesario tener las coordenadas y la cadena en la que se encuentra las regiones del genoma anotadas, para esto ensembl tiene estos archivos.       
-
-#se descarga el genoma de referencia con sus anotaciones, en este caso fueron obtenidos por ensembl100
-setwd("/home/aschafer/Documentos/Genomicas/Semestre_4/Genomica_humana/Enfermedades_monogenicas")
-homoS = import("Homo_sapiens.GRCh38.100.gff3.gz")
-head(homoS)
-
-#visualización de los datos que tiene el genoma de referencia. 
-table(mcols(homoS)$type)
-
-#nos quedamos solo con las columnas que deseamos.
-mcols(homoS) = mcols(homoS)[,c("source","type","ID","Name")]
-
-#Importamos nuestros datos a analizar
-Despues de haber trabajado nuestros datos con bowtie, utilizamos el archivo resultante .bam    
-
-#cargamos nuestros datos de las enfermedades monogenicas.
-bamFile <- "/home/aschafer/Documentos/Genomicas/Semestre_4/Genomica_humana/Enfermedades_monogenicas/sequences_aligned_sort.bam"
-informacion <- c("rname", "strand", "pos", "qwidth")
-informacion <- ScanBamParam(what = informacion)
-bam <- scanBam(bamFile, param = informacion)
-lapply(bam, names)
-
-#Construimos un GRanges
-mapGR <- GRanges(seqnames = bam[[1]]$rname, ranges = IRanges(start = bam[[1]]$pos, width = bam[[1]]$qwidth), strand = bam[[1]]$strand)
-head(mapGR)
-
-
-#Contamos mapeos dentro de las anotaciones 
-#buscamos en que regiones se traslapan 
-traslapados <- countOverlaps(homoS, mapGR)
-typeCounts <- aggregate(traslapados, by = list("type"=mcols(homoS)$type), sum)
-head(typeCounts)
-
-#vamos a graficar todos los resultados obtenidos, pero para eso vamos a juntas todos aquellos que tengan menos de 1000 juntos. 
-valMin <- 1000
-Cuentas <- typeCounts[typeCounts$x > valMin,]
-Cuentas <- Cuentas[order(Cuentas$x),]
-Cuentas_totales <- rbind(data.frame("type" = "other", "x" = sum(typeCounts$x[typeCounts$x <= valMin])), Cuentas)
-porcentaje <- round(100*Cuentas_totales$x/sum(Cuentas_totales$x), 1)
-
-pie(Cuentas_totales$x, labels = paste(porcentaje, "%", sep = ""), col = rev(rainbow(nrow(Cuentas_totales))), main="Porcentaje de lecturas totales alineadas por tipo", cex=0.9)
-legend("topleft", legend = Cuentas_totales$type, cex = 0.8, fill = rev(rainbow(nrow(Cuentas_totales))))
-
-Debido que vimos que los tipos: cromosoma, región biológica, exon, mRNA, gen, cinco prima UTR. No eran tan significativos, decidimos eliminarlos para los próximos análisis. 
-
-
-#graficamos ahora los nuevos datos, pero ahora en lugar de utilizar 1000 como el valor minimo, utilizamos 100.
-nuevosDatos <- typeCounts[typeCounts$type != "chromosome" & typeCounts$type != "biological_region" & typeCounts$type != "exon" & typeCounts$type != "mRNA" & typeCounts$type != "gene" & typeCounts$type != "five_prime_UTR",]
-colores <- c("violetred1", "dodgerblue4", "green2", "yellow", "red",  "steelblue1")
-valMin2 <- 100
-nuevas_cuentas <- nuevosDatos[nuevosDatos$x > valMin2,]
-nuevas_cuentas <- nuevas_cuentas[order(nuevas_cuentas$x),]
-nuevas_cuentas_totales <- rbind(data.frame("type"="other", "x" = sum(nuevosDatos$x[nuevosDatos$x <= valMin2])), nuevas_cuentas)
-porcentaje2 <- round(100*nuevas_cuentas_totales$x/sum(nuevas_cuentas_totales$x), 1)
-pie(nuevas_cuentas_totales$x, labels = paste(porcentaje2, "%", sep = ""), cex = 1, main = "Porcentaje de lecturas alineadas por tipo seleccionados", col = colores)
-legend("topleft", legend = nuevas_cuentas_totales$type, cex = 0.9, fill = colores)
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-