if(!require(keras)) {install.packages("keras"); library(keras)}
## The following method is also acceptable,
## But may cause it errors in proxy environments, etc.
## It uses R web settings to avoid Python web settings.
#
# mnist <- dataset_mnist()
#[1]
if(!file.exists("mnist.npz")){
  download.file("https://storage.googleapis.com/tensorflow/tf-keras-datasets/mnist.npz","mnist.npz")
}
mnist <- dataset_mnist(path=paste0(getwd(),"/mnist.npz"))
#preprocessing of data
x_train <- mnist$train$x; y_train <- mnist$train$y
x_test <- mnist$test$x; y_test <- mnist$test$y
x_train <- x_train / 255; x_test <- x_test / 255
## convert y matrix
num_classes <- 10
y_train <- to_categorical(y_train, num_classes)
y_test <- to_categorical(y_test, num_classes)
## convert x matrix
i_rows <- 28; i_cols <- 28
x_train <- array_reshape(x_train, c(nrow(x_train), i_rows, i_cols, 1))
x_test <- array_reshape(x_test, c(nrow(x_test), i_rows, i_cols, 1))
input_shape <- c(i_rows, i_cols, 1)
#[2]-----
MinistView <- function(data=x_test, range=1:16, col=4){
  par(mfrow=c(length(range)%/%col,col))
  par(mar=c(0.1,0.1,0.1,0.1))
  for(i in range){
    m <- matrix(data[i,,,1], 28,28)
    image(t(m[28:1,]), col=gray(255:0/255), axes=F)
  }
}
# It show x_miss 1:24 at 6 columns.
MinistView(x_train, 1:24, 6)
#[3]-----
# define model structure 
cnn_model <- keras_model_sequential() %>%
  layer_conv_2d(filters = 32, kernel_size = c(3,3), activation = 'relu', input_shape = input_shape) %>% 
  layer_max_pooling_2d(pool_size = c(2, 2)) %>% 
  layer_conv_2d(filters = 64, kernel_size = c(3,3), activation = 'relu') %>% 
  layer_max_pooling_2d(pool_size = c(2, 2)) %>% 
  layer_dropout(rate = 0.4) %>% 
  layer_flatten() %>% 
  layer_dense(units = 128, activation = 'relu') %>% 
  layer_dropout(rate = 0.4) %>% 
  layer_dense(units = num_classes, activation = 'softmax')
summary(cnn_model)
# Compile model
cnn_model %>% compile(
  loss = loss_categorical_crossentropy,
  optimizer = optimizer_adadelta(learning_rate = 1),
  metrics = c('accuracy')
)
# Train model
cnn_history <- cnn_model %>% fit(
  x_train, y_train,
  batch_size = 100,
  epochs = 15,
  validation_split = 0.2
)
cnn_model %>% evaluate(x_test, y_test)

## Check error images
cnn_pred <- cnn_model %>% predict(x_test) %>% k_argmax() %>% k_eval()
cat(paste("Number of misses: ",sum(cnn_pred != mnist$test$y),"/",length(y_test)))
y_miss <- mnist$test$y[cnn_pred != mnist$test$y]
pred_miss <- cnn_pred[cnn_pred != mnist$test$y]
x_miss <- mnist$test$x[cnn_pred != mnist$test$y,,]
x_miss <- array_reshape(x_miss, c(nrow(x_miss), i_rows, i_cols, 1))

# It show x_miss 1:16 at 16 columns.
MinistView(x_miss, 1:16, 16)
y_miss[1:16]
pred_miss[1:16]
