Runtime: 0.1375 sec
+Runtime: 0.137546 sec
Jackknife Statistics :
original bias std. error
- 100.029 100.019 0.150581
+ 99.9031 99.8931 0.149233
Polynomial degree: 1
+Polynomial degree: 1
Error: 0.08426840630693411
Bias^2: 0.07968918676726028
Var: 0.004579219539673833
0.08426840630693411 >= 0.07968918676726028 + 0.004579219539673833 = 0.08426840630693411
-Polynomial degree: 2
+Polynomial degree: 2
Error: 0.10398646080125035
Bias^2: 0.10077114273548986
Var: 0.0032153180657605086
@@ -868,21 +905,20 @@ Error: 0.05227921801205707
Bias^2: 0.048187277304303125
Var: 0.004091940707753964
0.05227921801205707 >= 0.048187277304303125 + 0.004091940707753964 = 0.05227921801205709
-Polynomial degree: 6
+Polynomial degree: 6
Error: 0.03781367141738898
Bias^2: 0.03365768507152761
Var: 0.004155986345861379
0.03781367141738898 >= 0.03365768507152761 + 0.004155986345861379 = 0.03781367141738899
-Polynomial degree: 7
+Polynomial degree:
+ 7
Error: 0.027609773491022498
Bias^2: 0.02299949826036597
Var: 0.004610275230656537
0.027609773491022498 >= 0.02299949826036597 + 0.004610275230656537 = 0.027609773491022505
-Polynomial degree: 8
+Polynomial degree: 8
Error: 0.017355848195591973
Bias^2: 0.010331721306655588
Var: 0.007024126888936384
@@ -892,7 +928,10 @@ Error: 0.026605727637189085
Bias^2: 0.010018312644140933
Var: 0.016587414993048166
0.026605727637189085 >= 0.010018312644140933 + 0.016587414993048166 = 0.0266057276371891
-Polynomial degree: 10
+Polynomial degree:
+ 10
Error: 0.021592704588043153
Bias^2: 0.010516485576652981
Var: 0.011076219011390184
@@ -902,23 +941,19 @@ Error: 0.07160048164228314
Bias^2: 0.01443680008897583
Var: 0.0571636815533073
0.07160048164228314 >= 0.01443680008897583 + 0.0571636815533073 = 0.07160048164228312
-Polynomial degree: 12
+Polynomial degree: 12
Error: 0.1154777721897675
Bias^2: 0.01628578269590588
Var: 0.09919198949386163
0.1154777721897675 >= 0.01628578269590588 + 0.09919198949386163 = 0.11547777218976751
-Polynomial degree: 13
+Polynomial degree: 13
Error: 0.22842468702166951
Bias^2: 0.01975416527163567
Var: 0.20867052175003387
0.22842468702166951 >= 0.01975416527163567 + 0.20867052175003387 = 0.22842468702166954
+
The bias-variance tradeoff summarizes the fundamental tension in @@ -1139,11 +1174,11 @@ Mean squared error on test data: 123711.53703498 Degree of polynomial: 3 Mean squared error on training data: 9011.85263220 Mean squared error on test data: 10913.84780262 -Degree of polynomial: 4 -Mean squared error on training data: 303.47610036 -
Mean squared error on test data: 426.30787294
+Degree of polynomial: 4
+Mean squared error on training data: 303.47610036
+Mean squared error on test data: 426.30787294
Degree of polynomial: 5
Mean squared error on training data: 3.80354994
Mean squared error on test data: 5.98822371
@@ -1175,64 +1210,62 @@ Mean squared error on test data: 0.17446471
Degree of polynomial: 13
Mean squared error on training data: 0.00759119
Mean squared error on test data: 1.08131003
-
-
-Degree of polynomial: 14
+Degree of polynomial: 14
Mean squared error on training data: 0.00472199
Mean squared error on test data: 0.81333793
-Degree of polynomial: 15
+
+
+Degree of polynomial: 15
Mean squared error on training data: 0.00410478
Mean squared error on test data: 92.09145189
Degree of polynomial: 16
Mean squared error on training data: 0.00315593
Mean squared error on test data: 234.39716546
-
-
-Degree of polynomial: 17
+Degree of polynomial: 17
Mean squared error on training data: 0.00242998
Mean squared error on test data: 1271.05295709
-Degree of polynomial: 18
+
+
+Degree of polynomial: 18
Mean squared error on training data: 0.00228740
Mean squared error on test data: 108.42208194
Degree of polynomial: 19
Mean squared error on training data: 0.00156372
Mean squared error on test data: 1388.41078073
-
-
-Degree of polynomial: 20
+Degree of polynomial: 20
Mean squared error on training data: 0.00137982
Mean squared error on test data: 1761.43341615
-Degree of polynomial: 21
+
+
+Degree of polynomial: 21
Mean squared error on training data: 0.00118170
Mean squared error on test data: 15061.31603087
Degree of polynomial: 22
Mean squared error on training data: 0.00092354
Mean squared error on test data: 890.63488525
-
-
-Degree of polynomial: 23
+Degree of polynomial: 23
Mean squared error on training data: 0.00085887
Mean squared error on test data: 5483.16796929
-Degree of polynomial: 24
+
+
+Degree of polynomial: 24
Mean squared error on training data: 0.00084589
Mean squared error on test data: 1695.57143061
Degree of polynomial: 25
Mean squared error on training data: 0.00078806
Mean squared error on test data: 131343.30655001
-
-
-Degree of polynomial: 26
+Degree of polynomial: 26
Mean squared error on training data: 0.00076916
Mean squared error on test data: 17709.14370264
-Degree of polynomial: 27
+
+
+Degree of polynomial: 27
Mean squared error on training data: 0.00068970
Mean squared error on test data: 2975.38903780
Degree of polynomial: 28
Mean squared error on training data: 0.00062588
Mean squared error on test data: 3848.64522721
-
-
-Degree of polynomial: 29
+Degree of polynomial: 29
Mean squared error on training data: 0.00060728
Mean squared error on test data: 2988.64001211
@@ -1243,7 +1276,7 @@ Mean squared error on test data: 2988.64001211
plt.plot(polynomial, np.log10(testerror), label='Test Error')
+
5.6. Exercises and Projects¶
+The main aim of this project is to study in more detail various +regression methods, including the Ordinary Least Squares (OLS) method, +The total score is 100 points. Each subtask has its own final score.
+We will first study how to fit polynomials to a specific +two-dimensional function called Franke’s +function. This +is a function which has been widely used when testing various +interpolation and fitting algorithms. Furthermore, after having +established the model and the method, we will employ resamling +techniques such as cross-validation and/or bootstrap in order to perform a +proper assessment of our models. We will also study in detail the +so-called Bias-Variance trade off.
+The Franke function, which is a weighted sum of four exponentials reads as follows
+The function will be defined for \(x,y\in [0,1]\). Our first step will +be to perform an OLS regression analysis of this function, trying out +a polynomial fit with an \(x\) and \(y\) dependence of the form \([x, y, +x^2, y^2, xy, \dots]\). We will also include bootstrap first as +a resampling technique. After that we will include the cross-validation technique. As in homeworks 1 and 2, we can use a uniform +distribution to set up the arrays of values for \(x\) and \(y\), or as in +the example below just a set of fixed +values for \(x\) and \(y\) with a given step +size. We will fit a +function (for example a polynomial) of \(x\) and \(y\). Thereafter we +will repeat much of the same procedure using the Ridge and Lasso +regression methods, introducing thus a dependence on the bias +(penalty) \(\lambda\).
+Finally we are going to use (real) digital terrain data and try to +reproduce these data using the same methods. We will also try to go +beyond the second-order polynomials metioned above and explore +which polynomial fits the data best.
+The Python code for the Franke function is included here (it performs also a three-dimensional plot of it)
+from mpl_toolkits.mplot3d import Axes3D
+import matplotlib.pyplot as plt
+from matplotlib import cm
+from matplotlib.ticker import LinearLocator, FormatStrFormatter
+import numpy as np
+from random import random, seed
+
+fig = plt.figure()
+ax = fig.gca(projection='3d')
+
+# Make data.
+x = np.arange(0, 1, 0.05)
+y = np.arange(0, 1, 0.05)
+x, y = np.meshgrid(x,y)
+
+
+def FrankeFunction(x,y):
+ term1 = 0.75*np.exp(-(0.25*(9*x-2)**2) - 0.25*((9*y-2)**2))
+ term2 = 0.75*np.exp(-((9*x+1)**2)/49.0 - 0.1*(9*y+1))
+ term3 = 0.5*np.exp(-(9*x-7)**2/4.0 - 0.25*((9*y-3)**2))
+ term4 = -0.2*np.exp(-(9*x-4)**2 - (9*y-7)**2)
+ return term1 + term2 + term3 + term4
+
+
+z = FrankeFunction(x, y)
+
+# Plot the surface.
+surf = ax.plot_surface(x, y, z, cmap=cm.coolwarm,
+ linewidth=0, antialiased=False)
+
+# Customize the z axis.
+ax.set_zlim(-0.10, 1.40)
+ax.zaxis.set_major_locator(LinearLocator(10))
+ax.zaxis.set_major_formatter(FormatStrFormatter('%.02f'))
+
+# Add a color bar which maps values to colors.
+fig.colorbar(surf, shrink=0.5, aspect=5)
+
+plt.show()
+
+5.6.1. Exercise: Ordinary Least Square (OLS) on the Franke function¶
+We will generate our own dataset for a function +\(\mathrm{FrankeFunction}(x,y)\) with \(x,y \in [0,1]\). The function +\(f(x,y)\) is the Franke function. You should explore also the addition +of an added stochastic noise to this function using the normal +distribution \(N(0,1)\).
+Write your own code (using either a matrix inversion or a singular +value decomposition from e.g., numpy ) or use your code from +homeworks 1 and 2 and perform a standard least square regression +analysis using polynomials in \(x\) and \(y\) up to fifth order. Find the +confidence intervals of the parameters (estimators) \(\beta\) by computing their +variances, evaluate the Mean Squared error (MSE)
+and the \(R^2\) score function. If \(\tilde{\hat{y}}_i\) is the predicted +value of the \(i-th\) sample and \(y_i\) is the corresponding true value, +then the score \(R^2\) is defined as
+where we have defined the mean value of \(\hat{y}\) as
+Your code has to include a scaling of the data (for example by +subtracting the mean value), and +a split of the data in training and test data. For this exercise you can +either write your own code or use for example the function for +splitting training data provided by the library Scikit-Learn (make +sure you have installed it). This function is called +\(train\_test\_split\). You should present a critical discussion of why and how you have scaled or not scaled the data.
+It is normal in essentially all Machine Learning studies to split the +data in a training set and a test set (eventually also an additional +validation set). There +is no explicit recipe for how much data should be included as training +data and say test data. An accepted rule of thumb is to use +approximately \(2/3\) to \(4/5\) of the data as training data.
+You can easily reuse the solutions to your exercises from week 35 and week 36.
+5.6.2. Exercise: Bias-variance trade-off and resampling techniques¶
+Our aim here is to study the bias-variance trade-off by implementing the bootstrap resampling technique.
+With a code which does OLS and includes resampling techniques, +we will now discuss the bias-variance trade-off in the context of +continuous predictions such as regression. However, many of the +intuitions and ideas discussed here also carry over to classification +tasks and basically all Machine Learning algorithms.
+Before you perform an analysis of the bias-variance trade-off on your test data, make +first a figure similar to Fig. 2.11 of Hastie, Tibshirani, and +Friedman. Figure 2.11 of this reference displays only the test and training MSEs. The test MSE can be used to +indicate possible regions of low/high bias and variance. You will most likely not get an +equally smooth curve!
+With this result we move on to the bias-variance trade-off analysis.
+Consider a +dataset \(\mathcal{L}\) consisting of the data +\(\mathbf{X}_\mathcal{L}=\{(y_j, \boldsymbol{x}_j), j=0\ldots n-1\}\).
+Let us assume that the true data is generated from a noisy model
+Here \(\epsilon\) is normally distributed with mean zero and standard +deviation \(\sigma^2\).
+In our derivation of the ordinary least squares method we defined then +an approximation to the function \(f\) in terms of the parameters +\(\boldsymbol{\beta}\) and the design matrix \(\boldsymbol{X}\) which embody our model, +that is \(\boldsymbol{\tilde{y}}=\boldsymbol{X}\boldsymbol{\beta}\).
+The parameters \(\boldsymbol{\beta}\) are in turn found by optimizing the means +squared error via the so-called cost function
+Here the expected value \(\mathbb{E}\) is the sample value.
+Show that you can rewrite this as
+Explain what the terms mean, which one is the bias and which one is +the variance and discuss their interpretations.
+Perform then a bias-variance analysis of the Franke function by +studying the MSE value as function of the complexity of your model.
+Discuss the bias and variance trade-off as function +of your model complexity (the degree of the polynomial) and the number +of data points, and possibly also your training and test data using the bootstrap resampling method.
+Note also that when you calculate the bias, in all applications you don’t know the function values \(f_i\). You would hence replace them with the actual data points \(y_i\).
+5.6.3. Exercise: Cross-validation as resampling techniques, adding more complexity¶
+The aim here is to write your own code for another widely popular +resampling technique, the so-called cross-validation method. Again, +before you start with cross-validation approach, you should scale your +data.
+Implement the \(k\)-fold cross-validation algorithm (write your own +code) and evaluate again the MSE function resulting +from the test folds. You can compare your own code with that from +Scikit-Learn if needed.
+Compare the MSE you get from your cross-validation code with the one +you got from your bootstrap code. Comment your results. Try \(5-10\) +folds. You can also compare your own cross-validation code with the +one provided by Scikit-Learn.
+5.6.4. Exercise: Ridge Regression on the Franke function with resampling¶
+Write your own code for the Ridge method, either using matrix +inversion or the singular value decomposition as done in the previous +exercise. Perform the same bootstrap analysis as in the +Exercise 2 (for the same polynomials) and the cross-validation in exercise 3 but now for different values of \(\lambda\). Compare and +analyze your results with those obtained in exercises 1-3. Study the +dependence on \(\lambda\).
+Study also the bias-variance trade-off as function of various values of +the parameter \(\lambda\). For the bias-variance trade-off, use the bootstrap resampling method. Comment your results.
+5.6.5. Exercise: Lasso Regression on the Franke function with resampling¶
+This exercise is essentially a repeat of the previous two ones, but now +with Lasso regression. Write either your own code (difficult and optional) or, in this case, +you can also use the functionalities of Scikit-Learn (recommended). +Give a +critical discussion of the three methods and a judgement of which +model fits the data best. Perform here as well an analysis of the bias-variance trade-off using the bootstrap resampling technique and an analysis of the mean squared error using cross-validation.
+5.6.6. Exercise: Analysis of real data¶
+With our codes functioning and having been tested properly on a +simpler function we are now ready to look at real data. We will +essentially repeat in this exercise what was done in exercises 1-5. However, we +need first to download the data and prepare properly the inputs to our +codes. We are going to download digital terrain data from the website +https://earthexplorer.usgs.gov/,
+Or, if you prefer, we have placed selected datafiles at https://github.com/CompPhysics/MachineLearning/tree/master/doc/Projects/2021/Project1/DataFiles
+In order to obtain data for a specific region, you need to register as +a user (free) at this website and then decide upon which area you want +to fetch the digital terrain data from. In order to be able to read +the data properly, you need to specify that the format should be SRTM +Arc-Second Global and download the data as a GeoTIF file. The +files are then stored in tif format which can be imported into a +Python program using
+scipy.misc.imread
+---------------------------------------------------------------------------
+NameError Traceback (most recent call last)
+<ipython-input-10-d985fb40c43d> in <module>
+----> 1 scipy.misc.imread
+
+NameError: name 'scipy' is not defined
+Here is a simple part of a Python code which reads and plots the data +from such files
+import numpy as np
+from imageio import imread
+import matplotlib.pyplot as plt
+from mpl_toolkits.mplot3d import Axes3D
+from matplotlib import cm
+
+# Load the terrain
+terrain1 = imread('SRTM_data_Norway_1.tif')
+# Show the terrain
+plt.figure()
+plt.title('Terrain over Norway 1')
+plt.imshow(terrain1, cmap='gray')
+plt.xlabel('X')
+plt.ylabel('Y')
+plt.show()
+If you should have problems in downloading the digital terrain data, +we provide two examples under the data folder of project 1. One is +from a region close to Stavanger in Norway and the other Møsvatn +Austfjell, again in Norway. +Feel free to produce your own terrain data.
+Alternatively, if you would like to use another data set, feel free to do so. This could be data close to your reseach area or simply a data set you found interesting. See for example kaggle.com for examples.
+Our final part deals with the parameterization of your digital terrain +data (or your own data). We will apply all three methods for linear regression, the same type (or higher order) of polynomial +approximation and cross-validation as resampling technique to evaluate which +model fits the data best.
+At the end, you should present a critical evaluation of your results +and discuss the applicability of these regression methods to the type +of data presented here (either the terrain data we propose or other data sets).
+