update exercises

doc/LectureNotes/exercisesweek37.ipynb

@@ -2,7 +2,7 @@
  "cells": [
   {
    "cell_type": "markdown",
-   "id": "ebfc874e",
+   "id": "cfebae03",
    "metadata": {
     "editable": true
    },
@@ -14,7 +14,7 @@
   },
   {
    "cell_type": "markdown",
-   "id": "8ecf95b0",
+   "id": "353ebb07",
    "metadata": {
     "editable": true
    },
@@ -27,7 +27,7 @@
   },
   {
    "cell_type": "markdown",
-   "id": "7940e437",
+   "id": "c02815dd",
    "metadata": {
     "editable": true
    },
@@ -45,7 +45,7 @@
   },
   {
    "cell_type": "markdown",
-   "id": "964cf4c4",
+   "id": "fb152043",
    "metadata": {
     "editable": true
    },
@@ -57,7 +57,7 @@
   },
   {
    "cell_type": "markdown",
-   "id": "bc8b093a",
+   "id": "9aba7a46",
    "metadata": {
     "editable": true
    },
@@ -68,7 +68,7 @@
   },
   {
    "cell_type": "markdown",
-   "id": "bc1907b5",
+   "id": "be52c097",
    "metadata": {
     "editable": true
    },
@@ -80,7 +80,7 @@
   },
   {
    "cell_type": "markdown",
-   "id": "c7d29897",
+   "id": "c54979e5",
    "metadata": {
     "editable": true
    },
@@ -90,7 +90,7 @@
   },
   {
    "cell_type": "markdown",
-   "id": "e45dbe57",
+   "id": "4b1eb811",
    "metadata": {
     "editable": true
    },
@@ -102,7 +102,7 @@
   },
   {
    "cell_type": "markdown",
-   "id": "0331d3dd",
+   "id": "036a5440",
    "metadata": {
     "editable": true
    },
@@ -114,7 +114,7 @@
   },
   {
    "cell_type": "markdown",
-   "id": "51159839",
+   "id": "f7fb8094",
    "metadata": {
     "editable": true
    },
@@ -125,7 +125,7 @@
   },
   {
    "cell_type": "markdown",
-   "id": "6df6ca33",
+   "id": "d7d92995",
    "metadata": {
     "editable": true
    },
@@ -137,7 +137,7 @@
   },
   {
    "cell_type": "markdown",
-   "id": "6bf812e3",
+   "id": "29eb9701",
    "metadata": {
     "editable": true
    },
@@ -150,7 +150,7 @@
   },
   {
    "cell_type": "markdown",
-   "id": "1e04714b",
+   "id": "aa106511",
    "metadata": {
     "editable": true
    },
@@ -162,17 +162,17 @@
   },
   {
    "cell_type": "markdown",
-   "id": "621afaa2",
+   "id": "58be1091",
    "metadata": {
     "editable": true
    },
    "source": [
-    "Show finally that the variance of $\\boldsymbol{\\beta}$ is"
+    "Show finally that the variance of $\\boldsymbol{\\boldsymbol{\\beta}}$ is"
    ]
   },
   {
    "cell_type": "markdown",
-   "id": "b7aa535d",
+   "id": "4cb58e41",
    "metadata": {
     "editable": true
    },
@@ -184,7 +184,7 @@
   },
   {
    "cell_type": "markdown",
-   "id": "45ea3a92",
+   "id": "a04df28d",
    "metadata": {
     "editable": true
    },
@@ -195,7 +195,7 @@
   },
   {
    "cell_type": "markdown",
-   "id": "3851ded2",
+   "id": "c095c96a",
    "metadata": {
     "editable": true
    },
@@ -207,7 +207,7 @@
   },
   {
    "cell_type": "markdown",
-   "id": "c759e617",
+   "id": "5d425a79",
    "metadata": {
     "editable": true
    },
@@ -220,7 +220,7 @@
   },
   {
    "cell_type": "markdown",
-   "id": "c16482df",
+   "id": "fb89cfab",
    "metadata": {
     "editable": true
    },
@@ -233,7 +233,7 @@
   },
   {
    "cell_type": "markdown",
-   "id": "2e2cda9b",
+   "id": "aea8801c",
    "metadata": {
     "editable": true
    },
@@ -245,7 +245,7 @@
   },
   {
    "cell_type": "markdown",
-   "id": "d69a45ab",
+   "id": "0e885693",
    "metadata": {
     "editable": true
    },

doc/src/week37/exercisesweek37.do.txt

@@ -53,7 +53,7 @@ With the OLS expressions for the optimal parameters $\bm{\hat{\beta}}$ show that
 \mathbb{E}(\bm{\hat{\beta}}) = \bm{\beta}.
 \]
 !et
-Show finally that the variance of $\bm{\beta}$ is
+Show finally that the variance of $\bm{\bm{\beta}}$ is
 !bt
 \[
 \mbox{Var}(\bm{\hat{\beta}}) = \sigma^2 \, (\mathbf{X}^{T} \mathbf{X})^{-1}.

doc/src/week37/programs/codeexamplesscaling.do.txt

@@ -5,7 +5,6 @@ DATE: today
 
 
 
-
 =====  This note contains code examples with a simple scaling  =====
 
 The programs here use both ordinrary least squares and Ridge regression with one value only for
@@ -40,12 +39,7 @@ true_beta = [2, 0.5, 3.7]
 
 # Make data set.
 x = np.linspace(-3, 3, n)
-y_real = 2 + 0.5*x + 3.7*x**2
-
-y = np.sum(
-    np.asarray([x ** p * b for p, b in enumerate(true_beta)]), 
-    axis=0) + 0.1 * np.random.normal(size=len(x))
-
+y = 2 + 0.5*x + 3.7*x**2
 
 #Design matrix X includes the intercept and scaling is made
 X = np.zeros((len(x), d))
@@ -216,22 +210,18 @@ Now we try to implement this.
 
 np.random.seed(2018)
 n = 100
-d = 3
+# we do not include the intercept
+d = 2
 Lambda = 0.01
 true_beta = [2, 0.5, 3.7]
 
 # Make data set.
 x = np.linspace(-3, 3, n)
-y_real = 2 + 0.5*x + 3.7*x**2
-
-y = np.sum(
-    np.asarray([x ** p * b for p, b in enumerate(true_beta)]), 
-    axis=0) + 0.1 * np.random.normal(size=len(x))
-
+y = 2 + 0.5*x + 3.7*x**2
 
 #Design matrix X does not include the intercept. 
 X = np.zeros((len(x), d))
-for p in range(d-1):     
+for p in range(d):     
     X[:, p] = x ** (p+1)
 
 
@@ -254,12 +244,13 @@ beta_OLS = OLS_fit_beta(X_train_scaled, y_train_scaled)
 beta_Ridge = Ridge_fit_beta(X_train_scaled, y_train_scaled,Lambda,d)
 print(beta_OLS)
 print(beta_Ridge)
-
+# calculate intercepts and print them
 interceptOLS = y_scaler - X_train_mean @ beta_OLS
 interceptRidge = y_scaler - X_train_mean @ beta_Ridge
 print(interceptOLS)
 print(interceptRidge)
-#predict value
+
+#predict value with intercept
 ytilde_test_OLS = X_test_scaled @ beta_OLS+y_scaler
 ytilde_test_Ridge = X_test_scaled @ beta_Ridge+y_scaler
 
@@ -275,7 +266,6 @@ print(MSE(y_test,ytilde_test_Ridge))
 
 
 plt.scatter(x,y,label='Data')
-#plt.plot(x,y_real,label='no noise')
 plt.plot(x, X @ beta_OLS+interceptOLS,'*', label="OLS_Fit")
 plt.plot(x, X @ beta_Ridge+interceptRidge, label="Ridge_Fit")
 plt.grid()
@@ -286,79 +276,70 @@ plt.show()
 
 We see that we get the same values for the parameters!  As it should be. The MSE may however change (not the case here).
 
-Finally, instead of using our own function we repeat the same example using the _standardscaler_ functionality of the library _Scikit-Learn_.
+Finally, instead of using our own function we repeat the same example
+using the _standardscaler_ functionality of the library
+_Scikit-Learn_.  Here we limit ourselves to Ridge regression only. 
 
-#!bc pycod
+!bc pycod
+
+from sklearn import linear_model
 np.random.seed(2018)
 n = 100
-d = 3
+d = 2
 Lambda = 0.01
 true_beta = [2, 0.5, 3.7]
 
 # Make data set.
 x = np.linspace(-3, 3, n)
-y_real = 2 + 0.5*x + 3.7*x**2
-
-y = np.sum(
-    np.asarray([x ** p * b for p, b in enumerate(true_beta)]), 
-    axis=0) + 0.1 * np.random.normal(size=len(x))
-
+y = (2 + 0.5*x + 3.7*x**2)
 
 #Design matrix X does not include the intercept. 
-X = np.zeros((len(x), d))
-for p in range(d-1):     
+X = np.zeros((n, d))
+for p in range(d):     
     X[:, p] = x ** (p+1)
 
-
 #Split data in train and test
 X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.2)
 
-# Scale data by subtracting mean value,own implementation
-#For our own implementation, we will need to deal with the intercept by centering the design matrix and the target variable
-X_train_mean = np.mean(X_train,axis=0)
-#Center by removing mean from each feature
-X_train_scaled = X_train - X_train_mean
-X_test_scaled = X_test - X_train_mean
-#The model intercept (called y_scaler) is given by the mean of the target variable (IF X is centered, note)
-y_scaler = np.mean(y_train)
-y_train_scaled = y_train - y_scaler
 
+# Scale data by subtracting mean value using scikit-learn
+from sklearn.preprocessing import StandardScaler
+scaler = StandardScaler()
+#scaler.fit(X_train)
+#scaler.fit(y_train)
+#X_train_scaled = scaler.transform(X_train)
+#X_test_scaled = scaler.transform(X_test)
+#y_train_scaled = scaler.transform(y_train)
 
 #Calculate beta
-beta_OLS = OLS_fit_beta(X_train_scaled, y_train_scaled)
-beta_Ridge = Ridge_fit_beta(X_train_scaled, y_train_scaled,Lambda,d)
-print(beta_OLS)
-print(beta_Ridge)
-
-interceptOLS = y_scaler - X_train_mean @ beta_OLS
-interceptRidge = y_scaler - X_train_mean @ beta_Ridge
-print(interceptOLS)
-print(interceptRidge)
-#predict value
-ytilde_test_OLS = X_test_scaled @ beta_OLS+y_scaler
-ytilde_test_Ridge = X_test_scaled @ beta_Ridge+y_scaler
-
+OLS = LinearRegression()
+OLS.fit(X_train,y_train)
+ypredictOLS = OLS.predict(X_test)
+RegRidge = linear_model.Ridge(Lambda)
+RegRidge.fit(X_train,y_train)
+ypredictRidge = RegRidge.predict(X_test)
+print(OLS.coef_)
+print(RegRidge.coef_)
+print(OLS.intercept_)
+interceptRidge = RegRidge.intercept_
+print(RegRidge.intercept_)
+#predict value without intercept
+ytilde_test_Ridge = X_test @ RegRidge.coef_+ RegRidge.intercept_
+ytilde_test_OLS = X_test @ OLS.coef_+ OLS.intercept_
 
 #Calculate MSE
-
 print("  ")
-print("test MSE of OLS:")
+print("test MSE of OLS")
 print(MSE(y_test,ytilde_test_OLS))
 print("  ")
 print("test MSE of Ridge")
 print(MSE(y_test,ytilde_test_Ridge))
-
-
 plt.scatter(x,y,label='Data')
-#plt.plot(x,y_real,label='no noise')
-plt.plot(x, X @ beta_OLS+interceptOLS,'*', label="OLS_Fit")
-plt.plot(x, X @ beta_Ridge+interceptRidge, label="Ridge_Fit")
+plt.plot(x, X @ RegRidge.coef_ + RegRidge.intercept_ , label="Ridge_Fit")
 plt.grid()
 plt.legend()
 plt.show()
-
-
-#!ec
+!ec