cnn.cpp (15035B)
1 #include "bpnn.hpp" 2 #include "utils.hpp" 3 4 //#include <Eigen/unsupported/CXX11/Tensor> 5 #include "cnn.hpp" 6 7 #define LARGE_NUM 1000000 // Remove me. 8 9 #if (!RECKLESS) 10 #define checknan(x, loc) if(x==INFINITY || x==NAN || x == -INFINITY) throw ValueError("Detected NaN in operation", loc) 11 #else 12 #define checknan(x, loc) 13 #endif 14 15 // NOTE: Below three functions not mine, from https://compvisionlab.wordpress.com/2014/01/01/c-code-for-reading-mnist-data-set/ 16 int ReverseInt (int i) 17 { 18 unsigned char ch1, ch2, ch3, ch4; 19 ch1=i&255; 20 ch2=(i>>8)&255; 21 ch3=(i>>16)&255; 22 ch4=(i>>24)&255; 23 return((int)ch1<<24)+((int)ch2<<16)+((int)ch3<<8)+ch4; 24 } 25 26 void ReadMNIST(int NumberOfImages, int DataOfAnImage,std::vector<std::vector<double>> &arr) 27 { 28 arr.resize(NumberOfImages,std::vector<double>(DataOfAnImage)); 29 std::ifstream file ("./t10k-images-idx3-ubyte",std::ios::binary); 30 if (file.is_open()) 31 { 32 int magic_number=0; 33 int number_of_images=0; 34 int n_rows=0; 35 int n_cols=0; 36 file.read((char*)&magic_number,sizeof(magic_number)); 37 magic_number= ReverseInt(magic_number); 38 file.read((char*)&number_of_images,sizeof(number_of_images)); 39 number_of_images= ReverseInt(number_of_images); 40 file.read((char*)&n_rows,sizeof(n_rows)); 41 n_rows= ReverseInt(n_rows); 42 file.read((char*)&n_cols,sizeof(n_cols)); 43 n_cols= ReverseInt(n_cols); 44 for(int i=0;i<number_of_images;++i) 45 { 46 for(int r=0;r<n_rows;++r) 47 { 48 for(int c=0;c<n_cols;++c) 49 { 50 unsigned char temp=0; 51 file.read((char*)&temp,sizeof(temp)); 52 arr[i][(n_rows*r)+c]= (double)temp; 53 } 54 } 55 } 56 } 57 } 58 59 unsigned char* read_mnist_labels(std::string full_path, int number_of_labels) { 60 auto reverseInt = [](int i) { 61 unsigned char c1, c2, c3, c4; 62 c1 = i & 255, c2 = (i >> 8) & 255, c3 = (i >> 16) & 255, c4 = (i >> 24) & 255; 63 return ((int)c1 << 24) + ((int)c2 << 16) + ((int)c3 << 8) + c4; 64 }; 65 66 typedef unsigned char uchar; 67 68 std::ifstream file(full_path, std::ios::binary); 69 70 if(file.is_open()) { 71 int magic_number = 0; 72 file.read((char *)&magic_number, sizeof(magic_number)); 73 magic_number = reverseInt(magic_number); 74 75 if(magic_number != 2049) throw std::runtime_error("Invalid MNIST label file!"); 76 77 file.read((char *)&number_of_labels, sizeof(number_of_labels)), number_of_labels = reverseInt(number_of_labels); 78 79 uchar* _dataset = new uchar[number_of_labels]; 80 for(int i = 0; i < number_of_labels; i++) { 81 file.read((char*)&_dataset[i], 1); 82 } 83 return _dataset; 84 } else { 85 throw std::runtime_error("Unable to open file `" + full_path + "`!"); 86 } 87 } 88 89 ConvLayer::ConvLayer(int x, int y, int stride, int kern_x, int kern_y, int pad, 90 std::function<float(float)> activ, std::function<float(float)> activ_deriv) 91 :stride_len(stride), padding(pad), activation(activ), activation_deriv(activ_deriv) 92 { 93 pad*=2; 94 input = new Eigen::MatrixXf (x+pad,y+pad); 95 dZ = new Eigen::MatrixXf (x+pad,y+pad); 96 for (int i = 0; i < (x+pad)*(y+pad); i++) { 97 (*input)((int)i / (y+pad),i%(y+pad)) = 0; 98 (*dZ)((int)i / (y+pad),i%(y+pad)) = 0; 99 } 100 kernel = new Eigen::MatrixXf (kern_x, kern_y); 101 for (int i = 0; i < kern_x*kern_y; i++) { 102 (*kernel)((int)i / kern_y,i%kern_y) = (float) rand() / RAND_MAX; 103 } 104 output = new Eigen::MatrixXf ((x-kern_x+1+pad/stride_len), (y-kern_y+1+pad/stride_len)); 105 for (int i = 0; i < (x-kern_y+1+pad/stride_len)*(y-kern_x+1+pad/stride_len); i++) { 106 (*output)((int)i / (y-kern_y+1+pad/stride_len),i%(y-kern_y+1+pad/stride_len)) = 0; 107 } 108 bias = 0; 109 }; 110 111 void ConvLayer::convolute() 112 { 113 for (int i = 0; i < (*input).rows(); i++) { 114 for (int j = 0; j < (*input).cols(); j++) { 115 (*dZ)(i, j) = activation_deriv((*input)(i, j)); 116 (*input)(i, j) = activation((*input)(i, j)); 117 } 118 } 119 for (int i = 0; i < output->rows(); i+=stride_len) { 120 for (int j = 0; j < output->cols(); j+=stride_len) { 121 (*output)(i, j) = (*kernel * (input->block(i, j, kernel->rows(), kernel->cols()))).sum(); 122 } 123 } 124 *output = (output->array() + bias).matrix(); 125 } 126 127 void ConvLayer::set_input(Eigen::MatrixXf* matrix) 128 { 129 input->block(padding, padding, matrix->rows(), matrix->cols()) = *matrix; 130 } 131 132 // Will eventually be different from ConvLayer 133 PoolingLayer::PoolingLayer(int x, int y, int stride, int kern_x, int kern_y, int pad) 134 :stride_len(stride), padding(pad) 135 { 136 input = new Eigen::MatrixXf (x+pad,y+pad); 137 for (int i = 0; i < (x+pad)*(y+pad); i++) { 138 (*input)((int)i / (y+pad),i%(y+pad)) = 0; 139 } 140 kernel = new Eigen::MatrixXf (kern_x, kern_y); 141 for (int i = 0; i < kern_x*kern_y; i++) { 142 (*kernel)((int)i / kern_y,i%kern_y) = (float) rand()/RAND_MAX; 143 } 144 output = new Eigen::MatrixXf (x-kern_x+1, y-kern_y+1); 145 for (int i = 0; i < (x-kern_x+1)*(y-kern_y+1); i++) { 146 (*output)((int)i / (y-kern_y+1),i%(y-kern_y+1)) = (float) rand()/RAND_MAX; 147 } 148 }; 149 150 void PoolingLayer::pool() 151 { 152 // It doesn't look like anything better than O(n^4) is doable for this as kernel needs to go through matrix and you need to index kernel. LOOK INTO ME!! 153 float maxnum = -LARGE_NUM; 154 for (int i = 0; i < input->cols() - kernel->cols(); i+=stride_len) { 155 for (int j = 0; j < input->rows() - kernel->rows(); j+=stride_len) { 156 for (int k = 0; k < kernel->cols(); k++) { 157 for (int l = 0; l < kernel->rows(); l++) { 158 if ((input->block(j, i, kernel->rows(), kernel->cols()))(l, k) > maxnum) { 159 maxnum = (input->block(j, i, kernel->rows(), kernel->cols()))(l, k); 160 } 161 } 162 } 163 } 164 } 165 } 166 167 ConvNet::ConvNet(const char* path, float learn_rate, float bias_rate, Regularization reg, float l, float ratio) 168 :Network(path, 1, learn_rate, bias_rate, reg, l, ratio), preprocess_length{0} 169 { 170 ReadMNIST(10000,784,data); 171 data_labels = read_mnist_labels("./t10k-labels-idx1-ubyte",10000); 172 labels = new Eigen::MatrixXf (1, 1); 173 } 174 175 void ConvNet::add_conv_layer(int x, int y, int stride, int kern_x, int kern_y, int pad, std::function<float(float)> activ, std::function<float(float)> activ_deriv) 176 { 177 preprocess_length+=1; 178 conv_layers.emplace_back(x,y,stride,kern_x, kern_y,pad,activ,activ_deriv); 179 } 180 181 // May make this inaccessible to user code and just have it called from add_conv_layer as pooling is basically always paired with conv. 182 void ConvNet::add_pool_layer(int x, int y, int stride, int kern_x, int kern_y, int pad) 183 { 184 pool_layers.emplace_back(x,y,stride,kern_x,kern_y,pad); 185 } 186 187 void ConvNet::initialize() 188 { 189 for (int i = 0; i < length-1; i++) { 190 layers[i].init_weights(layers[i+1]); 191 } 192 } 193 194 void ConvNet::next_batch() 195 { 196 for (int i = 0; i < 784; i++) { 197 (*conv_layers[0].input)(i/28, i%28) = data[batches][i]; 198 } 199 (*labels)(0,0) = (float)(int)data_labels[batches]; 200 } 201 202 void ConvNet::process() 203 { 204 // Assumes pooling is immediately after any conv layer. 205 for (int i = 0; i < preprocess_length-1; i++) { 206 conv_layers[i].convolute(); 207 // pool_layers[i].input = conv_layers[i].output; 208 // pool_layers[i].pool(); 209 conv_layers[i+1].input = conv_layers[i].output; 210 } 211 conv_layers[preprocess_length-1].convolute(); 212 //pool_layers[preprocess_length-1].input = conv_layers[preprocess_length-1].output; 213 //pool_layers[preprocess_length-1].pool(); 214 // std::cout << "Output:\n" << *pool_layers[preprocess_length-1].output << "\n\n"; 215 Eigen::Map<Eigen::RowVectorXf> flattened (conv_layers[preprocess_length-1].output->data(), 216 conv_layers[preprocess_length-1].output->size()); 217 // std::cout << "Flattened:\n" << flattened << "\n\n"; 218 for (int i = 0; i < flattened.cols(); i++) { 219 (*layers[0].contents)(0, i) = flattened[i]; 220 } 221 } 222 223 void ConvNet::set_label(Eigen::MatrixXf newlabels) 224 { 225 *labels = newlabels; 226 } 227 228 void ConvNet::list_net() 229 { 230 for (int i = 0; i < preprocess_length; i++) { 231 std::cout << "-----------------------\nCONVOLUTIONAL LAYER " << i 232 << "\n-----------------------\n\n\u001b[31mGENERAL INFO:\x1B[0;37m\nStride: " 233 << conv_layers[i].stride_len << "\nPadding: " << conv_layers[i].padding 234 << "\n\n\u001b[31mINPUT:\x1B[0;37m\n" << *conv_layers[i].input 235 << "\n\n\u001b[31mKERNEL:\x1B[0;37m\n" << *conv_layers[i].kernel 236 << "\n\n\u001b[31mOUTPUT:\x1B[0;37m\n" << *conv_layers[i].output 237 << "\n\n\u001b[31mBIAS:\x1B[0;37m\n" << conv_layers[i].bias << "\n\n\n"; 238 // std::cout << "-----------------------\nPOOLING LAYER " << i 239 // << "\n-----------------------\n\n\u001b[31mGENERAL INFO:\x1B[0;37m\nStride: " 240 // << pool_layers[i].stride_len << "\nPadding: " << conv_layers[i].padding 241 // << "\n\n\u001b[31mINPUT:\x1B[0;37m\n" << *pool_layers[i].input 242 // << "\n\n\u001b[31mKERNEL:\x1B[0;37m\n-" << *pool_layers[i].kernel 243 // << "\n\n\u001b[31mOUTPUT:\x1B[0;37m\n" << *pool_layers[i].output << "\n\n\n"; 244 } 245 std::cout << "-----------------------\nINPUT LAYER (LAYER 0)\n-----------------------" 246 << "\n\n\u001b[31mGENERAL INFO:\x1B[0;37m\nActivation Function: " 247 << layers[0].activation_str << "\n\n\u001b[31mACTIVATIONS:\x1B[0;37m\n" 248 << *layers[0].contents << "\n\n\u001b[31mWEIGHTS:\x1B[0;37m\n" 249 << *layers[0].weights << "\n\n\u001b[31mBIASES:\x1B[0;37m\n" 250 << *layers[0].bias << "\n\n\n"; 251 for (int i = 1; i < length-1; i++) { 252 std::cout << "-----------------------\nLAYER " << i 253 << "\n-----------------------\n\n\u001b[31mGENERAL INFO:\x1B[0;37m" 254 << "\nActivation Function: " << layers[i].activation_str 255 << "\n\n\u001b[31mACTIVATIONS:\x1B[0;37m\n" << *layers[i].contents 256 << "\n\n\u001b[31mBIASES:\x1B[0;37m\n" << *layers[i].bias 257 << "\n\n\u001b[31mWEIGHTS:\x1B[0;37m\n" << *layers[i].weights << "\n\n\n"; 258 } 259 std::cout << "-----------------------\nOUTPUT LAYER (LAYER " << length-1 260 << ")\n-----------------------\n\n\u001b[31mGENERAL INFO:\x1B[0;37m" 261 << "\nActivation Function: " << layers[length-1].activation_str 262 << "\n\n\u001b[31mACTIVATIONS:\x1B[0;37m\n" << *layers[length-1].contents 263 << "\n\n\u001b[31mBIASES:\x1B[0;37m\n" << *layers[length-1].bias << "\n\n\n"; 264 } 265 266 void ConvNet::backpropagate() 267 { 268 list_net(); 269 char a; 270 std::cin >> a; 271 std::vector<Eigen::MatrixXf> gradients; 272 gradients.push_back(Network::backpropagate()); 273 Eigen::Map<Eigen::MatrixXf> reshaped(gradients[gradients.size()-1].data(), 274 conv_layers.back().output->rows(), 275 conv_layers.back().output->cols()); 276 gradients[gradients.size()-1] = reshaped; 277 std::vector<Eigen::MatrixXf> conv_deltas; 278 for (int layer = conv_layers.size()-1; layer >= 0; layer--) { 279 conv_deltas.emplace_back(conv_layers[layer].kernel->rows(), 280 conv_layers[layer].kernel->cols()); 281 std::cout << conv_layers[layer].input->cols() << " " << gradients.back().cols() << "\n"; 282 for (int i = 0; i < conv_layers[layer].input->rows() - gradients.back().rows() + 1; i++) { 283 for (int j = 0; j < conv_layers[layer].input->cols() - gradients.back().cols() + 1; j++) { 284 conv_deltas[conv_deltas.size()-1](i, j) = 285 (gradients.back() * conv_layers[layer].input->block(i, j, gradients.back().rows(), 286 gradients.back().cols())).sum(); 287 } 288 } 289 *conv_layers[layer].kernel -= conv_deltas.back(); 290 291 Eigen::MatrixXf flipped_kernel = 292 Eigen::MatrixXf::Zero(conv_layers[layer].kernel->rows(), conv_layers[layer].kernel->cols()); 293 flipped_kernel = conv_layers[layer].kernel->transpose().colwise().reverse().transpose().colwise().reverse(); 294 Eigen::MatrixXf padded_grad = Eigen::MatrixXf::Zero(gradients.back().rows() + ((flipped_kernel.rows() - 1)*2), 295 gradients.back().cols() + ((flipped_kernel.cols() - 1)*2)); 296 padded_grad.block(flipped_kernel.rows() - 1, flipped_kernel.cols() - 1, 297 gradients.back().rows(), gradients.back().cols()) = gradients.back(); 298 Eigen::MatrixXf final_grad (padded_grad.rows() - flipped_kernel.rows() + 1, 299 padded_grad.cols() - flipped_kernel.cols() + 1); 300 for (int i = 0; i < padded_grad.rows() - flipped_kernel.rows() + 1; i++) { 301 for (int j = 0; j < padded_grad.cols() - gradients.back().cols() + 1; j++) { 302 final_grad(i, j) = (flipped_kernel * padded_grad.block(i, j, flipped_kernel.rows(), 303 flipped_kernel.cols())).sum(); 304 } 305 } 306 gradients.push_back(final_grad.cwiseProduct(*conv_layers[layer].dZ)); 307 } 308 list_net(); 309 assert(2<1); 310 } 311 312 void ConvNet::train() 313 { 314 float cost_sum = 0; 315 float acc_sum = 0; 316 for (int i = 0; i <= 100; i++) { 317 if (i != instances-batch_size) { // Don't try to advance batch on final batch. 318 next_batch(); 319 } 320 process(); 321 feedforward(); 322 backpropagate(); 323 cost_sum += cost(); 324 acc_sum += accuracy(); 325 batches++; 326 } 327 list_net(); 328 epoch_acc = 1.0/(100) * acc_sum; 329 epoch_cost = 1.0/(100) * cost_sum; 330 printf("Epoch %i complete - cost %f - acc %f\n", epochs, epoch_cost, epoch_acc); 331 batches=0; 332 decay(); 333 epochs++; 334 } 335 336 int main() 337 { 338 ConvNet net ("../data_banknote_authentication.txt", 0.05, 0.01, L2, 0, 0.9); 339 Eigen::MatrixXf labels (1,1); 340 net.add_conv_layer(28, 28, 1, 9, 9, 0, lecun_tanh, lecun_tanh_deriv); 341 // net.add_pool_layer(20,20,1,6,6,0); 342 net.add_conv_layer(20, 20, 1, 6, 6, 0, lecun_tanh, lecun_tanh_deriv); 343 std::cout << net.conv_layers[net.conv_layers.size()-1].output->rows() << "\n"; 344 //net.add_pool_layer(10,10,1,2,2,0); 345 net.add_layer(400, "sigmoid", sigmoid, sigmoid_deriv); 346 net.add_layer(5, "lecun_tanh", lecun_tanh, lecun_tanh_deriv); 347 net.add_layer(10, "resig", rectifier(sigmoid), rectifier(sigmoid_deriv)); 348 // net.list_net(); 349 // net.init_decay("step", 1, 2); 350 net.initialize(); 351 //net.list_net(); 352 353 for (int i = 0; i < 1; i++) { 354 net.train(); 355 } 356 net.list_net(); 357 }