R Notebook - Autoencoder Simples

Notebook com códigos R para a implementação de um Autoencoder simples - ERAMIA 2020

Pacote para trabalhar com imagens

library(png)

Carrega os digitos da base MNIST

carrega.digitos <- function(digitos = c(1,3,4,7,9), quantidade = 100,
                            treino=1, taxa.ruido = 0.1){

  digitos <<- digitos # Padroes a serem armazenados
  taxa.ruido <<- taxa.ruido # Probabilidade de um pixel ser trocado 
                              # (probabilidade de ruido) 

  ############## carrega os padroes #############
  padroes <<- list()
  orgdim <<- NULL

  cat("\nCarregando digitos...\n")

  for(i in digitos){

    if(treino == 1){

      arquivos <- list.files(paste("Digitos/MNIST/mnist_png/training/",
                                         as.character(i),"/",sep=""))[1:quantidade]

      for(j in 1:quantidade){
        img <- readPNG(paste("Digitos/MNIST/mnist_png/training/",
                                     as.character(i),"/",arquivos[j], sep=""))
        orgdim <<- dim(img)
        dim(img) <- NULL
        # Armazena os padroes
        padroes[[length(padroes)+1]] <<- as.vector(ifelse(img >= 0.2, 1, 0))
      }
    }
    else{
      arquivos <- list.files(paste("Digitos/MNIST/mnist_png/testing/",
                                   as.character(i),"/",sep=""))[1:quantidade]

      for(j in 1:quantidade){
        img <- readPNG(paste("Digitos/MNIST/mnist_png/testing/",
                                     as.character(i),"/",arquivos[j], sep=""))
        orgdim <<- dim(img)
        dim(img) <- NULL
        # Armazena os padroes
        padroes[[length(padroes)+1]] <<- as.vector(ifelse(img >= 0.2, 1, 0))
      }
    }
  }
}

Script para teste visual dos dígitos

testa.digitos <- function(modelo){

  # Carrega digitos de teste
  carrega.digitos(digitos = c(0,1,2,3,4,5,6,7,8,9),1,0,taxa.ruido)

  plotdim = 2*orgdim
  plot(c(1,(plotdim[1]+5)*length(digitos)), c(1,(plotdim[2]+5)*3),
            type="n", xlab="", ylab="")
  x = 1

  for (i in 1:length(padroes)) {
    padrao = padroes[[i]]

    ruido = (runif(length(padrao), 0, 1) > taxa.ruido ) * 1
    entrada = padrao * ruido
    retorno <- autoencoder.propagacao(modelo,entrada)
    ret <- retorno$y.escondida.saida
    ret <-  as.vector(ifelse(ret >= 0.5, 1, 0))

    # Padrao original
    img <- padrao;
    dim(img) <- orgdim
    image <- as.raster((img+1)/2)
    rasterImage(image, x, 1, x + plotdim[1], plotdim[2], interpolate=F)

    # Entrada com ruido
    img <- entrada;
    dim(img) <- orgdim
    image <- as.raster((img+1)/2)
    rasterImage(image, x, 1+(plotdim[2]+5), x + plotdim[1],
                1+2*(plotdim[2]+5), interpolate=F)

    # Imagem recuperada
    img <- ret;
    dim(img) <- orgdim
    image <- as.raster((img+1)/2)
    rasterImage(image, x, 1+2*(plotdim[2]+5), x + plotdim[1],
                1+2*(plotdim[2]+5)+plotdim[2], interpolate=F)

    x = x + plotdim[1]+5
  }
}

Vamos implementar uma função de ativação

funcao.ativacao <- function(v){

  # Função logística
  y <- 1 / (1 + exp(-v))

  return(y)
}

Vamos também precisar da derivada da função de ativação

der.funcao.ativacao <- function(y){

  # Derivada da logística
  derivada <- y * (1 - y)

  return(derivada)
}

Criação da arquitetura do Autoencoder. Similar à arquitetura do MLP, porém com uma camada de saída com a mesma quantidade de neurônios da camada de entrada

arquitetura <- function(num.entrada, num.escondida,
                        funcao.ativacao, der.funcao.ativacao){

  arq <- list()

  # Parametros da rede    
  arq$num.entrada.saida <- num.entrada
  arq$num.escondida <- num.escondida
  arq$funcao.ativacao <- funcao.ativacao
  arq$der.funcao.ativacao <- der.funcao.ativacao

  # 2 neuronios na camada escondida
  # 
  #       Ent1    Ent2   Bias 
  # 1     w11     w12    w13
  # 2     w21     w22    w23

  # Pesos conectando entrada e escondida
  num.pesos.entrada.escondida <- (num.entrada+1)*num.escondida
  arq$escondida <- matrix(runif(min=-0.5,max=0.5, num.pesos.entrada.escondida),
                                nrow=num.escondida, ncol=num.entrada+1)

  # Pesos conectando escondida e saida
  num.pesos.escondida.saida <- (num.escondida+1)*num.entrada
  arq$saida <- matrix(runif(min=-0.5,max=0.5,  num.pesos.escondida.saida),
                            nrow=num.entrada, ncol=num.escondida+1)

  return(arq)
}

Fase de propagação

autoencoder.propagacao <- function(arq, exemplo){

  # Entrada -> Cama Escondida
  v.entrada.escondida <- arq$escondida %*% as.numeric(c(exemplo,1))
  y.entrada.escondida <- arq$funcao.ativacao(v.entrada.escondida)

  # Camada Escondida -> Camada de Saida
  v.escondida.saida <- arq$saida %*% c(y.entrada.escondida,1)
  y.escondida.saida <- arq$funcao.ativacao(v.escondida.saida)

  # Resultados
  resultado <- list()
  resultado$v.entrada.escondida <- v.entrada.escondida
  resultado$y.entrada.escondida <- y.entrada.escondida
  resultado$v.escondida.saida <- v.escondida.saida
  resultado$y.escondida.saida <- y.escondida.saida

  return(resultado)
}

Fase de treinamento

autoencoder <- function(arq, dados, n, limiar){

  loss <- 10
  loss.anterior <- 0
  epocas <- 0

  # Treina eqto o erro for maior que um limiar
  while((abs(loss-loss.anterior)) > limiar){
        
    loss.anterior <- loss
    loss <- 0
        
    # Treino para todos os exemplos (epoca)
    for(i in 1:length(dados)){
            
      # Pego um exemplo de entrada
      x <- dados[[i]]
            
      # Pego a saida da rede para o exemplo
      resultado <- autoencoder.propagacao(arq,x)
      y <- resultado$y.escondida.saida
            
      # Calculo do erro para o exemplo
      erro <- x - y
      loss <- loss + -(sum(x*log(y) + (1-x)*log(1-y)))
            
      # Gradiente local no neuronio de saida
      # erro * derivada da funcao de ativacao
      grad.local.saida <- erro * arq$der.funcao.ativacao(y)
            
      # Gradiente local no neuronio escondido
      # derivada da funcao de ativacao no neuronio escondido * soma dos gradientes
      # locais dos neuronios conectados na proxima camada * pesos conectando a camada
            
      # escondida com a saida
      pesos.saida <- arq$saida[,1:arq$num.escondida]
      grad.local.escondida <- 
          as.numeric(arq$der.funcao.ativacao(resultado$y.entrada.escondida)) *
           (as.vector(grad.local.saida) %*% pesos.saida)
            
      # Ajuste dos pesos
      # Saida
      arq$saida <- arq$saida + n * (grad.local.saida %*% 
                                  c(resultado$y.entrada.escondida,1))
      # Escondida
      arq$escondida <- arq$escondida + n * (t(grad.local.escondida) %*%
                                                as.numeric(c(x,1)))
    }
  
    loss <- loss / length(dados)
    cat("Epoca = ",epocas," / Erro = ",loss," / ",abs(loss-loss.anterior),"\n")
    epocas <- epocas+1
  }

  retorno <- list()
  retorno$arq <- arq
  retorno$epocas <- epocas

  return(retorno)
}

Teste com camada subcompleta

# ===============================================================
# Carrega os digitos
carrega.digitos(digitos = c(0,1,2,3,4,5,6,7,8,9),50,1,0.3)

# Cria a arquitetuta
arq <- arquitetura(length(padroes[[1]]),length(padroes[[1]])-400,
                   funcao.ativacao,der.funcao.ativacao)

# Treina o modelo
modelo <- autoencoder(arq,padroes,0.2,0.1)

# Testa para novos digitos
testa.digitos(modelo$arq)
# ===============================================================

Teste com camada sobrecompleta

# ===============================================================
# Carrega os digitos
carrega.digitos(digitos = c(0,1,2,3,4,5,6,7,8,9),50,1,0.3)
#carrega.digitos2(digitos = c(0,1,2,3,4,5,6,7,8,9))

# Cria a arquitetuta
arq <- arquitetura(length(padroes[[1]]),length(padroes[[1]])+400,
                   funcao.ativacao,der.funcao.ativacao)

# Treina o modelo
modelo <- autoencoder(arq,padroes,0.2,0.1)

# Testa para novos digitos
testa.digitos(modelo$arq)
# ===============================================================