# # установка пакета mxnet. Достаточно откоментить и запустить только один раз, если mxnet ещё небыл установлен ранее.
# # Разработчики обещают выкладывать обновлённые версии каждую неделю, при желании можно обновлять пакет этими же командами.
#
# install.packages("drat", repos="https://cran.rstudio.com")
# drat:::addRepo("dmlc")
# install.packages("mxnet")
#
# # Есть версия пакета поддерживающая быстрые параллельные вычисления на видеокартах Nvidia.
# # Установка такой версии пакета сложнее, вот гайд - http://mxnet.io/get_started/setup.html#installing-mxnet-on-a-gpu
# # Используется технология CUDA, которая не поддерживается видюхами от AMD, очень плохо. 
# # При том что разработчики могли бы использовать OpenCL вместо CUDA, добившись быстроты и параллельности на всех видюхах, да ещё и на обычных CPU и плисах. Зачем они CUDA используют вообще хз, продались зелёным наверно.

library(mxnet)

# хелп по функциям:
# http://mxnet-mli.readthedocs.io/en/latest/packages/python/symbol.html#module-mxnet.symbol
# там хоть и для питона, но к R тоже много чего относится

# пример кода взят отсюдова: http://tjo-en.hatenablog.com/entry/2016/03/30/233848

#сколько бар будет использовано в каждом обучающем примере, и столько же будет использовано для предсказания
BARS_PER_SAMPLE <- 400
#количество обучающих примеров
TRAIN_SAMPLE_COUNT    <- 200
#количество тестовых примеров, они будут идти по времени после тренировочных
TEST_SAMPLE_COUNT     <- 200

ratesDF <- read.csv("D:\\eurusd_h1.csv", sep=";")

#считаем приросты цены по барам, для целевой переменной. diff уменьшает длинну вектора на 1, поэтому дописываем спереди NA чтоб вернуть первоначальную длинну
priceIncrease <- c(NA, diff(ratesDF[,"Open"]))

# цель переносим на 2 бара назад:
# 1) чтоб небыло заглядывания в будщее
# 2) потому что предсказываем используя бары 1-n а не бары 0-n (0 - текущий, несформированный, перерисовывающийся бар)
priceIncrease <- priceIncrease[c(-1,-2)]
#и убираем последних два бара из данных, а то на них и не обучиться, и ни протестировать
ratesDF <- ratesDF[c(-nrow(ratesDF), -(nrow(ratesDF)-1)), ]


#вектор с целевыми результатами
target <- priceIncrease
#три класса ("рост"(2), "без изменений"(1), "падение"(0)) вместо приростов цены.
target[target > 0] <- 2
target[target == 0] <- 1
target[target < 0] <- 0


#создаём матрицы для обучения и проверки модели

#mxnet использует вертикальные четырёхмерные таблицы для обучения
trainMatrix <- c()
for (endIndex in (nrow(ratesDF)-TRAIN_SAMPLE_COUNT-TEST_SAMPLE_COUNT+1):(nrow(ratesDF)-TEST_SAMPLE_COUNT)){
	startIndex <- endIndex-BARS_PER_SAMPLE+1
	
	#подавать просто цены в чистом виде нельзя, даже для свёрточной сети. Поэтому используем приращения цен.
	#Ещё на выбор можно нормализовывать каждый ряд цен в интервал (0:1) вместо использования приращений, но лично мне дельты кажутся более лучшими, стабильными. Кто знает опять таки, пробуйте всяко.
	
	#diff уменьшит длинну вектора на 1, так что нужен ряд цен на 1 длиннее чем итоговый вектор
	startIndex <- startIndex-1
	
	trainMatrix <- c(trainMatrix, diff(ratesDF[startIndex:endIndex,"Open"]))
	trainMatrix <- c(trainMatrix, diff(ratesDF[startIndex:endIndex,"High"]))
	trainMatrix <- c(trainMatrix, diff(ratesDF[startIndex:endIndex,"Low"]))
	trainMatrix <- c(trainMatrix, diff(ratesDF[startIndex:endIndex,"Close"]))
}
# конвертим вектор в 4-х мерную матрицу
# Измерения -
# 1) число бар на каждый обучающий пример
# 2) число векторов в каждом обучающем примере (тут их 4 - open, high, low, close)
# 3) параллельные каналы данных, если есть. Сейчас ничего нету, тут будет 1. Для обработки картинок были бы 3 канала - R, G, B
# 4) число обучающих примеров
# интуиция подсказывает что измерения 2 и 3 можно бы поменять местами, кажется логично, но это нужно ещё проверить. Для такого изменения нужно будет ещё и подправлять размеры кластеров в свёрточной сети.
dim(trainMatrix) <- c(BARS_PER_SAMPLE, 4, 1, TRAIN_SAMPLE_COUNT)
trainTarget <- target[(nrow(ratesDF)-TRAIN_SAMPLE_COUNT-TEST_SAMPLE_COUNT+1):(nrow(ratesDF)-TEST_SAMPLE_COUNT)]


testMatrix <- c()
for (endIndex in (nrow(ratesDF)-TEST_SAMPLE_COUNT+1):(nrow(ratesDF))){
	startIndex <- endIndex-BARS_PER_SAMPLE+1
	startIndex <- startIndex-1
	testMatrix <- c(testMatrix, diff(ratesDF[startIndex:endIndex,"Open"]))
	testMatrix <- c(testMatrix, diff(ratesDF[startIndex:endIndex,"High"]))
	testMatrix <- c(testMatrix, diff(ratesDF[startIndex:endIndex,"Low"]))
	testMatrix <- c(testMatrix, diff(ratesDF[startIndex:endIndex,"Close"]))
}
dim(testMatrix) <- c(BARS_PER_SAMPLE, 4, 1, TEST_SAMPLE_COUNT)
testTarget <- target[(nrow(ratesDF)-TEST_SAMPLE_COUNT+1):nrow(ratesDF)]






# создание структуры свёрточной сети
# 
# тут есть примеры структур для распознавания изображений, победители конкурсов и эстафет:
# https://adeshpande3.github.io/adeshpande3.github.io/The-9-Deep-Learning-Papers-You-Need-To-Know-About.html


#Входной слой
convNetStruct <- mx.symbol.Variable('data')

#первый свёрточный слой
#размер ядра 24x4 (24 бар, учитывая все 4 ценовых вектора OHLC), размер стоит подбирать разный, 24 бар это чисто наобум. Это просто офигеть какой важный параметр, и его надо подбирать
#Ещё не обязательно использовать все 4 вектора в ядре, можно попробовать с двумя или тремя.
#num_filter - число ядер. Опять таки никто не знает сколько надо, взято наобум.
convNetStruct <- mx.symbol.Convolution(data=convNetStruct, kernel=c(24,4), num_filter=round(BARS_PER_SAMPLE*0.5))
#активационная функция. На выбор есть "relu", "sigmoid", "softrelu", "tanh"
#relu вполне ок
#сигмоиду и тангенсоиду принято ругать и не использовать
#softrelu - пророчат светлое будущее, но таки кто знает... Стоит попробовать.
convNetStruct <- mx.symbol.Activation(data=convNetStruct, act_type="relu")
#группировка выходов по два, согласно максимального значения
convNetStruct <- mx.symbol.Pooling(data=convNetStruct, pool_type="max", kernel=c(2,1), stride=c(2,1))

# Второй свёрточный слой. Параметры и их смысл теже. Если ширина второго ядра была 4, то тут можно использовать максимум 1.
convNetStruct <- mx.symbol.Convolution(data=convNetStruct, kernel=c(48,1), num_filter=round(BARS_PER_SAMPLE*0.5^3))
convNetStruct <- mx.symbol.Activation(data=convNetStruct, act_type="relu")
convNetStruct <- mx.symbol.Pooling(data=convNetStruct, pool_type="max", kernel=c(2,1), stride=c(2,1))

# Таких свёрточных слоёв можно (нужно) делать десятки, иначе херня получится. Но для примера хватит и двух.

# переходим от свёрточных слоёв к обычным. Выходы последнего свёрточного слоя будут входами для обычной нейронки
convNetStruct <- mx.symbol.Flatten(data=convNetStruct)
# 30 нейронов в скрытом слое
convNetStruct <- mx.symbol.FullyConnected(data=convNetStruct, num_hidden=30)
# активационная функция
convNetStruct <- mx.symbol.Activation(data=convNetStruct, act_type="relu")

#выход
convNetStruct <- mx.symbol.FullyConnected(data=convNetStruct, num_hidden=3)
# Варианты активационной функции на выходе - LinearRegressionOutput, LogisticRegressionOutput, MAERegressionOutput, SoftmaxOutput
convNetStruct <- mx.symbol.SoftmaxOutput(data=convNetStruct)


#установка зерна ГПСЧ на 0, для получения одних и техже результатов каждый раз даже если перезапустить R. Это необязательно.
mx.set.seed(0)

# Тут нужно перекреститься и подождать часы/дни пока закончится обучение. Точность <= 0.5 для классификации означает полную бесполезность модели.
# У меня модель съела 4 gb оперативы чтоб обучиться даже с такой простой структурой, оперативки нужно очень много.
# num.round - число итераций обучения. С такими данными на обучение уйдут дни наверно. Поэтому тут просто 10 для непродолжительного обучения "чисто посмотреть".
#             Реально нужно поставить 10000 например, иногда приостанавливать обучение, смотреть ошибку на кроссвалидации, и полностью остановить обучение когда кроссвалидация начнёт стабильно ухудшаться.
# 
# общая инфа по параметрам нейронок - http://cs231n.github.io/neural-networks-3/
convModel <- mx.model.FeedForward.create(convNetStruct,
									X=trainMatrix, #предикторы
									y=trainTarget, #цель
									num.round=10,  #число итераций обучения
									array.batch.size = TRAIN_SAMPLE_COUNT, #число примеров используемых при коррекции весов на каждой итерации. Должно подойти количество самих обучающих примеров, но можно брать меньше для защиты от переобучения. Вроде будет хорошо если взять в два раза меньше чем число примеров
									optimizer = "sgd", #функция для изменения скорости обучения
									learning.rate = 0.0001, #начальная скорость обучения, она будет меняться в процессе обучения согласно параметру optimizer
									momentum = 0.5, #начальная инерция скорости обучения, она будет меняться в процессе обучения согласно параметру optimizer
									wd = 0.00001, #weight decay. Веса будут стремится к 0 с такой скоростью, чтоб избежать бесконечно увеличивающихся весов.
									eval.metric=mx.metric.accuracy, #модель стремится повысить точность
									epoch.end.callback=mx.callback.log.train.metric(100) #лог в консоль во время обучения
									)
# результат кстати плохой, 0.5. Нужна более сложная структура модели

# предсказание тренировочных и тестовых данных
trainPredictionMatrix <- predict(convModel, trainMatrix)
testPredictionMatrix <- predict(convModel, testMatrix)

trainPrediction <- c()
testPrediction  <- c()
for(i in 1:ncol(trainPredictionMatrix)){
	trainPrediction[i] <- which.max(trainPredictionMatrix[,i])-1
}
for(i in 1:ncol(testPredictionMatrix)){
	testPrediction[i] <- which.max(testPredictionMatrix[,i])-1
}

cat("Train classification accuracy: ", mean(trainPrediction == trainTarget), "\n")
cat("Test classification accuracy: ", mean(testPrediction == testTarget), "\n")