Esta practica traduce a R la misma logica que ya aprendiste en Python:
Todo el codigo usa solamente R base, para evitar problemas de paquetes.
2 + 2
## [1] 4
10 * 3
## [1] 30
25 / 5
## [1] 5
2^3
## [1] 8
En Python:
nombre = "TP53"
En R:
nombre <- "TP53"
organismo <- "Escherichia coli"
longitud <- 1500
gc <- 52.4
nombre
## [1] "TP53"
organismo
## [1] "Escherichia coli"
longitud
## [1] 1500
gc
## [1] 52.4
Podemos preguntar que tipo de dato contiene cada variable:
class(nombre)
## [1] "character"
class(longitud)
## [1] "numeric"
class(gc)
## [1] "numeric"
secuencia <- "ATGCGCGTAATCGGCTAA"
secuencia
## [1] "ATGCGCGTAATCGGCTAA"
nchar(secuencia)
## [1] 18
En Python usamos len(secuencia).
En R usamos nchar(secuencia) para contar caracteres.
bases <- strsplit(secuencia, split = "")[[1]]
bases
## [1] "A" "T" "G" "C" "G" "C" "G" "T" "A" "A" "T" "C" "G" "G" "C" "T" "A" "A"
conteo <- table(bases)
conteo
## bases
## A C G T
## 5 4 5 4
Podemos contar G y C de una forma muy segura:
g <- sum(bases == "G")
c <- sum(bases == "C")
g
## [1] 5
c
## [1] 4
g <- sum(bases == "G")
c <- sum(bases == "C")
longitud <- length(bases)
porcentaje_gc <- ((g + c) / longitud) * 100
porcentaje_gc
## [1] 50
round(porcentaje_gc, 2)
## [1] 50
secuencia_2 <- "GGGCGCGCATATATGCGCGC"
bases_2 <- strsplit(secuencia_2, split = "")[[1]]
g_2 <- sum(bases_2 == "G")
c_2 <- sum(bases_2 == "C")
gc_2 <- ((g_2 + c_2) / length(bases_2)) * 100
round(gc_2, 2)
## [1] 70
gc <- 58.8
if (gc > 50) {
print("GC alto")
} else {
print("GC bajo")
}
## [1] "GC alto"
Tres categorias:
gc <- 58.8
if (gc >= 60) {
categoria <- "GC alto"
} else if (gc >= 40) {
categoria <- "GC intermedio"
} else {
categoria <- "GC bajo"
}
categoria
## [1] "GC intermedio"
En Python trabajamos con listas.
En R una estructura fundamental es el vector.
valores_gc <- c(58.8, 55.3, 61.8)
muestras <- c("Integrante 1", "Integrante 2", "Integrante 3")
valores_gc
## [1] 58.8 55.3 61.8
muestras
## [1] "Integrante 1" "Integrante 2" "Integrante 3"
En R el primer elemento tiene indice 1:
valores_gc[1]
## [1] 58.8
Estadistica basica:
mean(valores_gc)
## [1] 58.63333
median(valores_gc)
## [1] 58.8
min(valores_gc)
## [1] 55.3
max(valores_gc)
## [1] 61.8
datos <- data.frame(
muestra = c("Integrante 1", "Integrante 2", "Integrante 3"),
organismo = c(
"Aurantimonas litoralis",
"Pseudomonas putida",
"Catenulispora acidiphila"
),
gc = c(58.8, 55.3, 61.8)
)
datos
## muestra organismo gc
## 1 Integrante 1 Aurantimonas litoralis 58.8
## 2 Integrante 2 Pseudomonas putida 55.3
## 3 Integrante 3 Catenulispora acidiphila 61.8
Explorar la tabla:
head(datos)
## muestra organismo gc
## 1 Integrante 1 Aurantimonas litoralis 58.8
## 2 Integrante 2 Pseudomonas putida 55.3
## 3 Integrante 3 Catenulispora acidiphila 61.8
str(datos)
## 'data.frame': 3 obs. of 3 variables:
## $ muestra : chr "Integrante 1" "Integrante 2" "Integrante 3"
## $ organismo: chr "Aurantimonas litoralis" "Pseudomonas putida" "Catenulispora acidiphila"
## $ gc : num 58.8 55.3 61.8
Seleccionar una columna:
datos$gc
## [1] 58.8 55.3 61.8
Promedio de la columna:
mean(datos$gc)
## [1] 58.63333
Filtrar:
datos[datos$gc > 58, ]
## muestra organismo gc
## 1 Integrante 1 Aurantimonas litoralis 58.8
## 3 Integrante 3 Catenulispora acidiphila 61.8
for (valor in valores_gc) {
print(valor)
}
## [1] 58.8
## [1] 55.3
## [1] 61.8
Con una condicion:
for (valor in valores_gc) {
if (valor >= 60) {
mensaje <- "GC alto"
} else {
mensaje <- "GC menor de 60"
}
print(paste(valor, mensaje))
}
## [1] "58.8 GC menor de 60"
## [1] "55.3 GC menor de 60"
## [1] "61.8 GC alto"
En Python construimos funciones con def.
En R:
calcular_gc <- function(seq) {
bases <- strsplit(seq, split = "")[[1]]
g <- sum(bases == "G")
c <- sum(bases == "C")
gc <- ((g + c) / length(bases)) * 100
return(gc)
}
Probamos nuestra funcion:
calcular_gc("ATGCGCGTAATCGGCTAA")
## [1] 50
round(calcular_gc("ATGCGCGTAATCGGCTAA"), 2)
## [1] 50
secuencias <- c(
"ATGCGCGTAATCGGCTAA",
"ATATATATATGCGCGCAA",
"GGGCGCGCGCGCATATAT"
)
secuencias
## [1] "ATGCGCGTAATCGGCTAA" "ATATATATATGCGCGCAA" "GGGCGCGCGCGCATATAT"
Aplicamos la misma funcion a todas:
gc_resultados <- sapply(secuencias, calcular_gc)
gc_resultados
## ATGCGCGTAATCGGCTAA ATATATATATGCGCGCAA GGGCGCGCGCGCATATAT
## 50.00000 33.33333 66.66667
Creamos una tabla:
resultados <- data.frame(
muestra = c("Muestra_1", "Muestra_2", "Muestra_3"),
secuencia = secuencias,
gc = round(gc_resultados, 2)
)
resultados
## muestra secuencia gc
## ATGCGCGTAATCGGCTAA Muestra_1 ATGCGCGTAATCGGCTAA 50.00
## ATATATATATGCGCGCAA Muestra_2 ATATATATATGCGCGCAA 33.33
## GGGCGCGCGCGCATATAT Muestra_3 GGGCGCGCGCGCATATAT 66.67
barplot(
resultados$gc,
names.arg = resultados$muestra,
ylab = "Porcentaje GC",
main = "Contenido GC por muestra"
)
Agregar el promedio:
barplot(
resultados$gc,
names.arg = resultados$muestra,
ylab = "Porcentaje GC",
main = "Contenido GC por muestra"
)
abline(
h = mean(resultados$gc),
lty = 2
)
Ahora usamos una tabla sencilla para conectar R basico con transcriptomica.
expresion <- data.frame(
gen = c("BRCA1", "TP53", "MYC", "GAPDH"),
control = c(12.4, 8.1, 5.2, 40.0),
tratamiento = c(31.7, 7.9, 18.6, 39.5)
)
expresion
## gen control tratamiento
## 1 BRCA1 12.4 31.7
## 2 TP53 8.1 7.9
## 3 MYC 5.2 18.6
## 4 GAPDH 40.0 39.5
Calculamos el cambio:
expresion$cambio <- expresion$tratamiento - expresion$control
expresion
## gen control tratamiento cambio
## 1 BRCA1 12.4 31.7 19.3
## 2 TP53 8.1 7.9 -0.2
## 3 MYC 5.2 18.6 13.4
## 4 GAPDH 40.0 39.5 -0.5
Gen con mayor aumento:
expresion[which.max(expresion$cambio), ]
## gen control tratamiento cambio
## 1 BRCA1 12.4 31.7 19.3
Grafico:
matriz_exp <- as.matrix(
expresion[, c("control", "tratamiento")]
)
barplot(
t(matriz_exp),
beside = TRUE,
names.arg = expresion$gen,
legend.text = c("Control", "Tratamiento"),
ylab = "Expresion",
main = "Control vs tratamiento"
)
Esto NO es todavia un analisis de RNA-seq.
Es una introduccion conceptual al tipo de tablas que despues
analizaremos con herramientas como DESeq2.
secuencias_final <- c(
M1 = "ATGCGCGTAATCGGCTAA",
M2 = "ATATATATATATATATGC",
M3 = "GGGCGCGCGCGCATATAT",
M4 = "ATGCCCGGGTTTAAACCC"
)
gc_final <- sapply(secuencias_final, calcular_gc)
gc_final
## M1 M2 M3 M4
## 50.00000 11.11111 66.66667 55.55556
resultados_final <- data.frame(
muestra = names(secuencias_final),
gc = round(gc_final, 2)
)
resultados_final
## muestra gc
## M1 M1 50.00
## M2 M2 11.11
## M3 M3 66.67
## M4 M4 55.56
resultados_final[which.max(resultados_final$gc), ]
## muestra gc
## M3 M3 66.67
mean(resultados_final$gc)
## [1] 45.835
barplot(
resultados_final$gc,
names.arg = resultados_final$muestra,
ylab = "Porcentaje GC",
main = "Contenido GC de las secuencias"
)
| Concepto | Python | R |
|---|---|---|
| Variable | x = 5 |
x <- 5 |
| Texto | "ATGC" |
"ATGC" |
| Longitud de texto | len(x) |
nchar(x) |
| Lista / vector | [1,2,3] |
c(1,2,3) |
| Primer elemento | x[0] |
x[1] |
| Condicion | if x > 5: |
if (x > 5) {} |
| Ciclo | for x in lista: |
for (x in vector) {} |
| Funcion | def funcion(x): |
funcion <- function(x) {} |
| Tabla | DataFrame |
data.frame |
| Promedio | mean() |
mean() |
| Grafico | matplotlib |
barplot() |
Si entiendes estas tres ideas:
datos$gc
## [1] 58.8 55.3 61.8
datos[datos$gc > 58, ]
## muestra organismo gc
## 1 Integrante 1 Aurantimonas litoralis 58.8
## 3 Integrante 3 Catenulispora acidiphila 61.8
mean(datos$gc)
## [1] 58.63333
calcular_gc <- function(seq) {
# instrucciones
return(resultado)
}
sapply(secuencias, calcular_gc)
ya reconoces tres pilares de R:
trabajar con tablas + crear funciones + automatizar analisis.
La sintaxis cambia respecto a Python.
La logica que ya aprendiste permanece.