diff --git a/MLPP/LinReg/LinReg.cpp b/MLPP/LinReg/LinReg.cpp index faf2c1e..8a42857 100644 --- a/MLPP/LinReg/LinReg.cpp +++ b/MLPP/LinReg/LinReg.cpp @@ -34,6 +34,36 @@ namespace MLPP{ return Evaluate(x); } + void LinReg::NewtonRaphson(double learning_rate, int max_epoch, bool UI){ + LinAlg alg; + Reg regularization; + double cost_prev = 0; + int epoch = 1; + forwardPass(); + while(true){ + cost_prev = Cost(y_hat, outputSet); + + std::vector error = alg.subtraction(y_hat, outputSet); + + // Calculating the weight gradients (2nd derivative) + std::vector first_derivative = alg.mat_vec_mult(alg.transpose(inputSet), error); + std::vector> second_derivative = alg.matmult(alg.transpose(inputSet), inputSet); + weights = alg.subtraction(weights, alg.scalarMultiply(learning_rate/n, alg.mat_vec_mult(alg.transpose(alg.inverse(second_derivative)), first_derivative))); + weights = regularization.regWeights(weights, lambda, alpha, reg); + + // Calculating the bias gradients (2nd derivative) + bias -= learning_rate * alg.sum_elements(error) / n; // We keep this the same. The 2nd derivative is just [1]. + forwardPass(); + + if(UI) { + Utilities::CostInfo(epoch, cost_prev, Cost(y_hat, outputSet)); + Utilities::UI(weights, bias); + } + epoch++; + if(epoch > max_epoch) { break; } + } + } + void LinReg::gradientDescent(double learning_rate, int max_epoch, bool UI){ LinAlg alg; Reg regularization; @@ -59,7 +89,6 @@ namespace MLPP{ Utilities::UI(weights, bias); } epoch++; - if(epoch > max_epoch) { break; } } } diff --git a/MLPP/LinReg/LinReg.hpp b/MLPP/LinReg/LinReg.hpp index cdb58f8..56b5ef8 100644 --- a/MLPP/LinReg/LinReg.hpp +++ b/MLPP/LinReg/LinReg.hpp @@ -17,6 +17,7 @@ namespace MLPP{ LinReg(std::vector> inputSet, std::vector outputSet, std::string reg = "None", double lambda = 0.5, double alpha = 0.5); std::vector modelSetTest(std::vector> X); double modelTest(std::vector x); + void NewtonRaphson(double learning_rate, int max_epoch, bool UI); void gradientDescent(double learning_rate, int max_epoch, bool UI = 1); void SGD(double learning_rate, int max_epoch, bool UI = 1); void MBGD(double learning_rate, int max_epoch, int mini_batch_size, bool UI = 1); diff --git a/MLPP/NumericalAnalysis/NumericalAnalysis.cpp b/MLPP/NumericalAnalysis/NumericalAnalysis.cpp index 52f6332..3c348a0 100644 --- a/MLPP/NumericalAnalysis/NumericalAnalysis.cpp +++ b/MLPP/NumericalAnalysis/NumericalAnalysis.cpp @@ -5,6 +5,7 @@ // #include "NumericalAnalysis.hpp" +#include "LinAlg/LinAlg.hpp" #include namespace MLPP{ @@ -20,6 +21,18 @@ namespace MLPP{ return (function(x + eps) - 2 * function(x) + function(x - eps)) / (eps * eps); } + double NumericalAnalysis::constantApproximation(double(*function)(double), double c){ + return function(c); + } + + double NumericalAnalysis::linearApproximation(double(*function)(double), double c, double x){ + return constantApproximation(function, c) + numDiff(function, c) * (x - c); + } + + double NumericalAnalysis::quadraticApproximation(double(*function)(double), double c, double x){ + return linearApproximation(function, c, x) + 0.5 * numDiff_2(function, c) * (x - c) * (x - c); + } + double NumericalAnalysis::numDiff(double(*function)(std::vector), std::vector x, int axis){ // For multivariable function analysis. // This will be used for calculating Jacobian vectors. @@ -77,4 +90,18 @@ namespace MLPP{ } return hessian; } + + double NumericalAnalysis::constantApproximation(double(*function)(std::vector), std::vector c){ + return function(c); + } + + double NumericalAnalysis::linearApproximation(double(*function)(std::vector), std::vector c, std::vector x){ + LinAlg alg; + return constantApproximation(function, c) + alg.matmult(alg.transpose({jacobian(function, c)}), {alg.subtraction(x, c)})[0][0]; + } + + double NumericalAnalysis::quadraticApproximation(double(*function)(std::vector), std::vector c, std::vector x){ + LinAlg alg; + return linearApproximation(function, c, x) + 0.5 * alg.matmult({(alg.subtraction(x, c))}, alg.matmult(hessian(function, c), alg.transpose({alg.subtraction(x, c)})))[0][0]; + } } \ No newline at end of file diff --git a/MLPP/NumericalAnalysis/NumericalAnalysis.hpp b/MLPP/NumericalAnalysis/NumericalAnalysis.hpp index ee3313f..d0f5c2b 100644 --- a/MLPP/NumericalAnalysis/NumericalAnalysis.hpp +++ b/MLPP/NumericalAnalysis/NumericalAnalysis.hpp @@ -18,6 +18,10 @@ namespace MLPP{ double numDiff(double(*function)(double), double x); double numDiff_2(double(*function)(double), double x); + double constantApproximation(double(*function)(double), double c); + double linearApproximation(double(*function)(double), double c, double x); + double quadraticApproximation(double(*function)(double), double c, double x); + double numDiff(double(*function)(std::vector), std::vector x, int axis); double numDiff_2(double(*function)(std::vector), std::vector x, int axis1, int axis2); @@ -25,6 +29,11 @@ namespace MLPP{ std::vector jacobian(double(*function)(std::vector), std::vector x); // Indeed, for functions with scalar outputs the Jacobians will be vectors. std::vector> hessian(double(*function)(std::vector), std::vector x); + + double constantApproximation(double(*function)(std::vector), std::vector c); + double linearApproximation(double(*function)(std::vector), std::vector c, std::vector x); + double quadraticApproximation(double(*function)(std::vector), std::vector c, std::vector x); + }; } diff --git a/README.md b/README.md index 59b1c53..f1a7e9e 100644 --- a/README.md +++ b/README.md @@ -149,7 +149,11 @@ The result will be the model's predictions for the entire dataset. - Multivariate Functions 2. Jacobian Vector Calculator 3. Hessian Matrix Calculator - 3. Newton-Raphson Method + 4. Function approximator + - Constant Approximation + - Linear Approximation + - Quadratic Approximation + 5. Newton-Raphson Method 14. ***Linear Algebra Module*** 15. ***Statistics Module*** 16. ***Data Processing Module*** diff --git a/main.cpp b/main.cpp index a847ee0..6bb65f5 100644 --- a/main.cpp +++ b/main.cpp @@ -49,8 +49,12 @@ using namespace MLPP; +// double f(double x){ +// return x*x*x + 2*x - 2; +// } + double f(double x){ - return x*x*x + 2*x - 2; + return cos(x); } /* y = x^3 + 2x - 2 @@ -77,7 +81,7 @@ double f_mv(std::vector x){ Where x, y = x[0], x[1], this function is defined as: f(x, y) = x^4 + xy^3 + yz^2 - ∂f/∂x = 3x^2 + 3y^2 + ∂f/∂x = 4x^3 + 3y^2 ∂^2f/∂x^2 = 6x ∂f/∂z = 2zy @@ -178,17 +182,37 @@ int main() { // UniLinReg model(inputSet, outputSet); // alg.printVector(model.modelSetTest(inputSet)); - // MULIVARIATE LINEAR REGRESSION + // // MULIVARIATE LINEAR REGRESSION // std::vector> inputSet = {{1,2,3,4,5,6,7,8,9,10}, {3,5,9,12,15,18,21,24,27,30}}; // std::vector outputSet = {2,4,6,8,10,12,14,16,18,20}; + // LinReg model(alg.transpose(inputSet), outputSet); // Can use Lasso, Ridge, ElasticNet Reg - // model.normalEquation(); + // model.gradientDescent(0.001, 30000, 1); // model.SGD(0.001, 30000, 1); // model.MBGD(0.001, 10000, 2, 1); + // model.normalEquation(); + // alg.printVector(model.modelSetTest((alg.transpose(inputSet)))); // std::cout << "ACCURACY: " << 100 * model.score() << "%" << std::endl; + + // std::cout << "Total epoch num: 300" << std::endl; + // std::cout << "Method: 1st Order w/ Jacobians" << std::endl; + + // LinReg model(alg.transpose(inputSet), outputSet); // Can use Lasso, Ridge, ElasticNet Reg + + // model.gradientDescent(0.001, 300, 0); + + + // std::cout << "--------------------------------------------" << std::endl; + // std::cout << "Total epoch num: 300" << std::endl; + // std::cout << "Method: Newtonian 2nd Order w/ Hessians" << std::endl; + // LinReg model2(alg.transpose(inputSet), outputSet); + + // model2.NewtonRaphson(1.5, 300, 0); + + // // LOGISTIC REGRESSION // std::vector> inputSet; // std::vector outputSet; @@ -474,8 +498,6 @@ int main() { // data.getImage("../../Data/apple.jpeg", chicken); // alg.printVector(chicken); - // // TESTING QR DECOMP. EXAMPLE VIA WIKIPEDIA. SEE https://en.wikipedia.org/wiki/QR_decomposition. - // std::vector> P = {{12, -51, 4}, {6, 167, -68}, {-4, 24, -41}}; // alg.printMatrix(P); @@ -504,6 +526,10 @@ int main() { // Checks for numerical analysis class. NumericalAnalysis numAn; + //std::cout << numAn.quadraticApproximation(f, 0, 1) << std::endl; + + std::cout << numAn.quadraticApproximation(f_mv, {0, 0, 0}, {1, 1, 1}) << std::endl; + // std::cout << numAn.numDiff(&f, 1) << std::endl; // std::cout << numAn.newtonRaphsonMethod(&f, 1, 1000) << std::endl; @@ -513,9 +539,9 @@ int main() { //std::cout << numAn.numDiff_2(&f, 2) << std::endl; - std::cout << numAn.numDiff_2(&f_mv, {2, 2, 500}, 2, 2) << std::endl; - std::cout << "Our Hessian." << std::endl; - alg.printMatrix(numAn.hessian(&f_mv, {2, 2, 500})); + // std::cout << numAn.numDiff_2(&f_mv, {2, 2, 500}, 2, 2) << std::endl; + // std::cout << "Our Hessian." << std::endl; + // alg.printMatrix(numAn.hessian(&f_mv, {2, 2, 500})); return 0; }