Cómo Programar una Red Neuronal desde 0 en R

Share on linkedin
Share on twitter
Share on email

En R existen varios paquetes que nos permiten crear Redes Neuronales, como neuralnet o las más recientes (y conocidas) tensorflow y keras. Sí son paquetes muy potentes y seguramente escriba sobre ellos en un futuro, pero me gustaría comenzar explicando qué es una red neuronal, qué partes tiene y cómo podemos programar una desde 0 en R. Así que vamos a ello.

Datos para nuestra Red Neuronal

Lo primero que necesitamos es un problema a resolver. En nuestro caso, siguiendo el ejemeplo de DotCSV, crearemos un problema de clasificación.

Para ello, he creado una función que genera datos circulares según un radio R. Con esto, podremos dibujar dos círculos concéntricos que se solapen un poco, para ver qué tal se le da a nuestra red neuronal.

Esto lo he hecho con el siguiente código y siguiendo este enlace.

circulo <- function(x, R, centroX=0, centroY=0){
r = R * sqrt(runif(x))
theta = runif(x) * 2 * pi
x = centroX + r * cos(theta)
y = centroY + r * sin(theta)

z = data.frame(x = x, y = y)
return(z)
}
datos1 <- circulo(150,0.5)
datos2 <- circulo(150,1.5)

datos1$Y <- 1
datos2$Y <- 0
datos <- rbind(datos1,datos2)

rm(datos1,datos2, circulo)

Así pues, nuestros datos tienen la siguiente forma:

library(ggplot2)
ggplot(datos,aes(x,y, col = as.factor(Y))) + geom_point()

X <- as.matrix(datos[,1:2])
Y <- as.matrix(datos[,3])

rm(datos)

Creamos las matrices X e Y y con esto ya tenemos nuestro problema de clasificación listo para ser resuelto con nuestra red neuronal programada desde cero en R. ¡Vamos a ello!

Partes de una red neuronal

La idea es sencilla: programar cada parte de la red por separado, para así ir explicando en qué consisten las redes neuronales. Una vez tengamos todo programado, probaremos la red neuronal en nuestro problema para ver que funciona correctamente.

Así pues, podríamos decir que una red neuronal se programa en 3 partes:

  • La estructura de la red neuronal: simplemente defines el número de capas, neuronas por cada capa, funciones de activación, la capa de salida… en definitiva, el esqueleto. Además, damos unos valores iniciales a dichos parámetros.

  • Forward Propagation: “arrancamos” la red neuronal y hacemos que nos devuelva un resultado.

  • Back Propagation: en base al error de la estimación, ajustamos los parámetros de la última capa a las capas anteriores. Además, usamos el descenso del gradiente (Gradient Descent) para reajustar los parámetros e ir entrenando a la red neuronal.

Por último tendríamos que entrenar a la red en un problema real y ver cómo funciona. Fácil, ¿verdad? Pues al lío 😉

Programar una red neuronal en R: la estructura de la red neuronal

Las redes neuronales son neuronas conectadas entre ellas de manera secuencial. A cada una de esas secuencias las llamamos capas.

Las capas siempre tienen las misma estructura: neuronas de la capa, número de conexiones, parámetros b (bias) y w… Por eso, como siempre es la misma estruco código más simple.

La forma en la que lo he hecho ha sido mediante objetos S4. No hay mucha info al respecto… pero bueno, ( en este enlace puedes aprender cómo funcionan).

Esta clase tendrá los siguientes tres parámetros:

  • Número de neuronas en la capa.
  • Número de neuronas en la capa anterior.
  • Bias(b): el número que se suma a la neurona antes de a pasarla por la función de activación.
  • El número de conexiones que entran a la capa (W): indica el número de conexiones que entran en la capa. Como cada neurona está conectada con todas las neuronas de la capa siguiente, el número de conexiones de una capa dependerá del número de neuronas en la capa anterior y el número de neuronas en la capa actual.
  • Función de activación dentro de la capa: indicando la “transformación” que se aplicará a en las neuronas de esa capa. En este caso, todas las neuronas dentro de la misma capa aplicarán la misma función de activación, aunque se podrán aplicar una función diferente en cada capa.

Así pues, creamos una clase S4 que contenga estos elementos. Además, aprovecharemos esta clase para inicializar los valores W y b de forma aleatoria, usando la función runif.

neurona <- setRefClass(
  "neurona",
  fields = list(
    fun_act = "list",
    numero_conexiones = "numeric",
    numero_neuronas = "numeric",
    W = "matrix",
    b = "numeric"
  ),
  methods = list(
    initialize = function(fun_act, numero_conexiones, numero_neuronas)
    {
      fun_act <<- fun_act
      numero_conexiones <<- numero_conexiones
      numero_neuronas <<- numero_neuronas
      W <<- matrix(runif(numero_conexiones*numero_neuronas),
                   nrow = numero_conexiones)
      b <<- runif(numero_neuronas)
    }
  )
)

Asimismo, deberemos programar las funciones de activación para poder usarlas en la red neuronal. En este caso, usaremos dos tipos de funciones de activación: función sigmoide y función relu.

Además, aprovecharemos estas funciones, no solo para que nos devuelvan el valor de la función en sí, sino también su derivada. ¿Por qué? Pues porque las derivadas son necesarias para calcular el descenso del gradiente.

Además, incluir los dos valores en la misma función es más sencillo a la hora de programar. ¿Que quieres el resultado de la función de activación? Eliges el primer valor de la lista que devuelve. ¿Que quieres la derivada? Eliges el segundo valor.

Así pues, empezamos programac con la función sigmoide, que es así:

sigmoid = function(x) {
  y = list() 
  y[[1]] <- 1 / (1 + exp(-x))
  y[[2]] <- x * (1 - x)
  return(y)
}

x <- seq(-5, 5, 0.01)
plot(x, sigmoid(x)[[1]], col='blue')

Ahora solo nos queda programar la función relu, junto con su derivada, como hago a continuación:

relu <- function(x){
  y <- list()
  y[[1]] <- ifelse(x<0,0,x)
  y[[2]] <- ifelse(x<0,0,1)
  return(y)
}

plot(x, relu(x)[[1]], col='blue')

Por último, nos queda poner todo de forma conjunta para crear la estructura de la red neuronal, es decir, las capas que conforman la red.

En nuestro caso, al haber creado la clase neurona, crearemos las capas de manera iterativa. Así, esta misma forma de programar las neuronas nos servirá también para otras redes neurionales con otra estructura.

n = ncol(X) #Núm de neuronas en la primera capa
capas = c(n, 4, 8, 1) # Número de neuronas en cada capa.
funciones = list(sigmoid, relu, sigmoid) # Función de activación en cada capa

red <- list()

for (i in 1:(length(capas)-1)){
    red[[i]] <- neurona$new(funciones[i],capas[i], capas[i+1])
}

red
## [[1]]
## Reference class object of class "neurona"
## Field "fun_act":
## [[1]]
## function(x) {
##   y = list() 
##   y[[1]] <- 1 / (1 + exp(-x))
##   y[[2]] <- x * (1 - x)
##   return(y)
## }
## 
## Field "numero_conexiones":
## [1] 2
## Field "numero_neuronas":
## [1] 4
## Field "W":
##           [,1]       [,2]      [,3]       [,4]
## [1,] 0.8927513 0.01044008 0.1007646 0.02930838
## [2,] 0.9394562 0.90888047 0.5629580 0.05452293
## Field "b":
## [1] 0.001307318 0.491861450 0.043629507 0.646223775
## 
## [[2]]
## Reference class object of class "neurona"
## Field "fun_act":
## [[1]]
## function(x){
##   y <- list()
##   y[[1]] <- ifelse(x<0,0,x)
##   y[[2]] <- ifelse(x<0,0,1)
##   return(y)
## }
## 
## Field "numero_conexiones":
## [1] 4
## Field "numero_neuronas":
## [1] 8
## Field "W":
##           [,1]       [,2]       [,3]      [,4]      [,5]      [,6]      [,7]
## [1,] 0.1344001 0.93911520 0.16928257 0.3929345 0.4333714 0.9781074 0.5214181
## [2,] 0.2581032 0.89714491 0.06853257 0.7858252 0.1900454 0.6251393 0.7785054
## [3,] 0.7606003 0.03149109 0.29002985 0.1246016 0.2214262 0.5777163 0.9180025
## [4,] 0.9574416 0.23177728 0.89906514 0.6796380 0.2312405 0.8389011 0.7654668
##           [,8]
## [1,] 0.2547775
## [2,] 0.1115950
## [3,] 0.9771444
## [4,] 0.5807271
## Field "b":
## [1] 0.01464125 0.02918942 0.92919122 0.02876864 0.39818770 0.10067225 0.89313200
## [8] 0.92939015
## 
## [[3]]
## Reference class object of class "neurona"
## Field "fun_act":
## [[1]]
## function(x) {
##   y = list() 
##   y[[1]] <- 1 / (1 + exp(-x))
##   y[[2]] <- x * (1 - x)
##   return(y)
## }
## 
## Field "numero_conexiones":
## [1] 8
## Field "numero_neuronas":
## [1] 1
## Field "W":
##            [,1]
## [1,] 0.97304955
## [2,] 0.05584097
## [3,] 0.62309299
## [4,] 0.27844900
## [5,] 0.45473066
## [6,] 0.67233558
## [7,] 0.86945649
## [8,] 0.30373311
## Field "b":
## [1] 0.08716024

Ahora que ya tenemos la estructura de la red neuronal, pero de momento esto nos sirve de poco. Tenemos que entranar a nuestra red neuronal, así que vamos a ello.

Programar una red neuronal en R: entrenando a la red neuronal

El entrenamiento de una red neuronal tiene tres pasos:

  • Front propagation: la red neuronal recibirá unos valores de entrada. Capa a capa se irán haciendo unas operaciones hasta que la red nos devuelva un valor. Seguramente, como la red no esté entrenada fallará más que una escopeta de feria, por lo que habrá que ir haciendo ajustes, mediante los dos siguientes pasos:

  • Back propagation: una vez hemos recibido el resultado, lo comparamos con el resultado real. Esta comparación la hacemos aplicando una función de coste, lo cual nos devolverá un error. Utilizaremos ese error para ir ajustando los valores de la red neuronal hacia atrás: última capa, anteúltima, etc.

  • Optimizamos la función de coste: a medida que se va propagando hacia atrás, usamos una función de optimización (generalemente el descenso del gradiente) para ir optimizando los parámetros b y W. Así, la red irá mejorando en el entrenamiento poco a poco.

Hechas las introducciones, ¡vamos con la propagación hacia adelante!

Entrenando a la Red Neuronal: Propagación hacia adelante

A la hora de hacer la propagación hacia adelante, iremos guardando el resultado de la neurona antes (a) y después de aplicar la función de activación (z). ¿Por qué? Pues para tener los valores guardados para cuando hagamos la propagación hacia atrás, más que nada.

entrenar <- function(red, X,Y, coste){

out = list()
out[[1]] <- append(list(matrix(0,ncol=2,nrow=1)), list(X))
  
  for(i in c(1:(length(red)))){
    z = list((out[[length(out)]][[2]] %*% red[[i]]$W + red[[i]]$b))
    a = list(red[[i]]$fun_act[[1]](z[[1]])[[1]])
    out[[i+1]] <- append(z,a)
  }
return(out)
}
#rm(a,out,z,i)

Si ejectuamos la función, tendremos los valores predichos, como vemos a continuación:

forward <- entrenar(red, X,Y, coste)
head(forward[[4]][[2]])
##           [,1]
## [1,] 0.9949141
## [2,] 0.9992319
## [3,] 0.9995250
## [4,] 0.9992261
## [5,] 0.9984179
## [6,] 0.9994190

Pero claro, esta predicción… pues seguramente esté lejos de ser certera. Y es que, tenemos que entrenar a nuestra red neuronal. Pero, ¿cómo sabemos cuánto falla? Para eso existe la función de coste;)

Entrenando a la Red Neuronal: Función de Coste

La función de coste es una función que compara todas las predicciones con el valor real para saber en qué medida nuestra red neuronal ha fallado y así poder ir mejorándola.

En nuestro caso vamos a usar una función de coste muy típica: el error cuadrático medio, que es el resultado de hacer la media de las diferencias entre el valor predicho y el valor real elevados al cuadrado.

Esta es una función de coste típica que se usa en muchos algortimos, como la regresión lineal. Básicamente la idea es que lo que falla el algoritmo es la distancia al cuadrado entre los puntos predichos sobre el valor real.

¿Por qué la distancia al cuadrado? Pues para que los errores positivos (valor predicho > valor real, es decir, te has “pasado” en la predicción) y los errores negativos (te has quedado “corto” en la predicción) no se neutralicen entre ellos.

Al igual que hemos hecho con las funciones de activación, aprovecharemos la para que la función de coste nos devuelva su derivada parcial, puesto que la necesitaremos para el descenso del gradiente.

coste <- function(Yp,Yr){
  y <- list()
  y[[1]] <- mean((Yp-Yr)^2)
  y[[2]] <- (Yp-Yr)
  return(y)
}

Ya tenemos nuestra función de coste para “juzgar” nuestra red neuronal programada en R. Ahora, vamos a entrenarla.

Entrenando a la Red Neuronal: Back Propagation

En este caso, tenemos que hacer la propagación hacia atrás. Para ello, usaremos la funciónrev, la cual devuelve el reverso del objeto que le pasemos, en nuestro caso, la lista out, que son los resultados por capa del front propagation.

Ahora, en cada caso tendremos que calcular el error de la capa respecto al coste, también conocido como la delta. Para calcular la delta respecto a la capa previa, habría dos casos:

  • En la última capa debemos multiplicar la derivada de la función de coste y el resultado de la la función de activación en la última capa.

  • En las demás capas, debemos usar el error (delta) de la capa anterior, multiplicado de forma matricial por el vector de pesos de la capa siguiente, para así moverlo hacia atrás.

delta <- list() 

for (i in rev(1:length(red))){
  z = out[[i+1]][[1]]
  a = out[[i+1]][[2]]
  
  if(i == length(red)){
    delta[[1]] <- coste(a,Y)[[2]] * red[[i]]$fun_act[[1]](a)[[2]]
  } else{
    delta <- list(delta[[1]] %*% t(red[[i]]$W) * red[[i]]$fun_act[[1]](a)[[2]],delta)
  }
}

¡Backpropagation listo! Vamos ahora con el punto final: gradient descent.

Entrenando a la Red Neuronal: Gradient Descent

Por último, vamos a optimizar los parámetros de nuestra red mediante las deltas que ya tenemos calculadas. Como gradient descent requiere pasar por cada capa para optimizar los parámetros, vamos a aprovechar el pase hacia atrás que ya hemos programado.

Así pues, lo que haremos será actualizar los valores de by W. Además, usaremos un ratio conocido como learning Rate para determinar en cómo queremos que la red entrene.

Sin embargo, el learning rate es algo sensible ya que si es muy pequeño, el aprendizaje será más “seguro”, pero mucho más lento, mientras que si es muy alto, puede que la función de activación no converja. Os dejo un enlace a un vídeo de Deeplearning.ai donde Andrew Ng lo explica de forma muy clara.

En cualquier caso, la función de backpropagation quedaría así:

backprop <- function(out, red, lr = 0.05){

  delta <- list() 
  
  for (i in rev(1:length(red))){
    z = out[[i+1]][[1]]
    a = out[[i+1]][[2]]
    
    if(i == length(red)){
      delta[[1]] <- coste(a,Y)[[2]] * red[[i]]$fun_act[[1]](a)[[2]]
    } else{
      delta <- list(delta[[1]] %*% W_temp * red[[i]]$fun_act[[1]](a)[[2]],delta)
    }
    
    W_temp = t(red[[i]]$W)
    
    red[[i]]$b <- red[[i]]$b - mean(delta[[1]]) * lr
    red[[i]]$W <- red[[i]]$W - t(out[[i]][[2]]) %*% delta[[1]] * lr
  }
  x = list()
  x[[1]] <- red
  x[[2]] <- out[[length(out)]][[1]]
  return(x)
  
}

Juntando todo: entrenando a la red neuronal desde 0

Ya hemos programado la red neuronal desde 0 en R. Ahora queda lo mejor: ver cómo funciona.

Para ello, lo primero que vamos a hacer es juntar la parte de front propagation con la de backpropagation y gradient descent en una misma función:

red_neuronal <- function(red, X,Y, coste,lr = 0.05){
## Front Prop
out = list()
out[[1]] <- append(list(matrix(0,ncol=4,nrow=1)), list(X))
  
  for(i in c(1:(length(red)))){
    z = list((out[[length(out)]][[2]] %*% red[[i]]$W + red[[i]]$b))
    a = list(red[[i]]$fun_act[[1]](z[[1]])[[1]])
    out[[i+1]] <- append(z,a)
  }

  
## Backprop & Gradient Descent
delta <- list() 
  
for (i in rev(1:length(red))){
  z = out[[i+1]][[1]]
  a = out[[i+1]][[2]]
    
  if(i == length(red)){
    delta[[1]] <- coste(a,Y)[[2]] * red[[i]]$fun_act[[1]](a)[[2]]
  } else{
    delta <- list(delta[[1]] %*% W_temp * red[[i]]$fun_act[[1]](a)[[2]],delta)
  }
    
  W_temp = t(red[[i]]$W)
    
  red[[i]]$b <- red[[i]]$b - mean(delta[[1]]) * lr
  red[[i]]$W <- red[[i]]$W - t(out[[i]][[2]]) %*% delta[[1]] * lr
  
}
return(out[[length(out)]][[2]])

}

Vamos a probar que nuestra función red neuronal funciona correctamente.

resultado <- red_neuronal(red, X,Y, coste)
dim(resultado)
## [1] 300   1

¡Funciona! Ahora solo queda lo último. Iterar para que la red vaya aprendiendo. Además, iremos mostrando el error de la red para ir viendo como va mejorando a medida que pasa el tiempo.

Para ello, vamos a hacer 1000 iteraciones e ir calculando y guardando el error cada 20 iteraciones.

for(i in seq(1000)){

  Yt = red_neuronal(red, X,Y, coste, lr=0.01)

  if(i %% 25 == 0){
    if(i == 25){
      iteracion <- i
      error <- coste(Yt,Y)[[1]]
    }else{
          iteracion <- c(iteracion,i)
          error <- c(error,coste(Yt,Y)[[1]])      
    }
  }
}

Y por último, visualizamos el error de nuestra red neuronal:

library(ggplot2)

grafico = data.frame(Iteracion = iteracion,Error = error)

ggplot(grafico,aes(iteracion, error)) + geom_line() + theme_minimal() +
  labs(title = "Evolución del error de la Red Neuronal")

¡Funciona! Como véis, la red neuronal ha ido aprendiendo poco a poco hasta llegar a las 750 iteraciones, donde ya ha dejado de aprender, porque la red ya está optimizada.

Asi que esto es todo. Espero te haya servido y que hayas podido aprender nuevas cosas. ¡Nos vemos en el próximo post!