Problem Set 2 — Advanced Data Analytics

HEC Lausanne, autumn 2026 · Topics: pandas, multiple linear regression, Newton's method

Based on course materials by Prof. Simon Scheidegger.

How to work through this notebook. Cells marked TODO are yours to complete — code cells contain ... where your code goes, markdown cells say (your answer here). Everything runs top to bottom with no internet access and no absolute paths; the data ships in data/. When you are done, Kernel → Restart & Run All must execute cleanly to the end.

Nothing is submitted. These problem sets are practice; the solution is discussed in the exercise session. Because nobody grades them, the notebook has to tell you whether you got it right — that is what the "Certify it" parts are for: checks you run on your own code, with the expected outcome stated so you can see for yourself.

AI assistants. You may — and realistically will — work with an AI assistant (ChatGPT, Claude, Copilot). That is fine. Assistants are very good at producing code that runs and looks right; the skill that will carry your project is telling whether it is right. Treat every "Certify it" part as the question you should be asking the assistant's code — and your own.

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from sklearn.linear_model import LinearRegression

Exercise 1 — Multiple regression: CO₂ emissions of cars

Multiple regression is linear regression with more than one explanatory variable. The dataset data/cars.csv lists 36 car models with their engine size (Volume, in cm³), weight (Weight, in kg), and CO₂ emissions (CO2, in g/km).

We want to predict a car's CO₂ emissions from its weight and engine volume:

$$\mathrm{CO_2} = \beta_0 + \beta_1\,\mathrm{Weight} + \beta_2\,\mathrm{Volume} + \varepsilon$$

(a) Fit and predict

Load the data, build the feature matrix X (columns Weight, Volume) and the target y (CO2), and fit a sklearn.linear_model.LinearRegression. Then predict the CO₂ emission of a car that weighs 2300 kg and has a 1300 cm³ engine, and print the result.

df = pd.read_csv("data/cars.csv")

# TODO: build X (columns Weight, Volume) and y (column CO2),
# fit a LinearRegression, and predict CO2 for Weight=2300, Volume=1300.
# Keep X as a DataFrame so the model remembers the feature names.
...

(b) Coefficients

Print the fitted coefficients and the intercept. Then answer in the markdown cell below: by how much does predicted CO₂ change when the weight increases by 1 kg, holding engine volume fixed?

# TODO: print the coefficients (with their feature names) and the intercept.
...

TODO (your answer here)

(c) Certify it

The model is linear, so it must obey a property you can check exactly: the effect of +1000 kg is $1000 \cdot \beta_1$, no matter which car you start from. Contract check: predict CO₂ for a car of 3300 kg and 1300 cm³ and assert that the difference from your part-(a) prediction equals $1000 \cdot \beta_1$ to within 1e-9. If your assistant wrote part (a) for you, this is where you find out whether it fitted what you think it fitted.

# TODO: predict for Weight=3300, Volume=1300; print the new prediction and
# assert that (new - old) == 1000 * beta_Weight (to within 1e-9).
...

Exercise 2 — Newton's method for a nonlinear system

Solve the following system with a self-implemented Newton's method, as derived in the lecture — no scipy.optimize for the iteration itself:

$$\begin{aligned} f_1(x_1, x_2) &= 5x_1^2 - x_2^2 = 0\\ f_2(x_1, x_2) &= x_2 - \tfrac{1}{4}\left(\sin x_1 + \cos x_2\right) = 0 \end{aligned}$$

Newton's method for systems: given an iterate $x^{(k)}$, solve the linear system $J(x^{(k)})\, \delta = -f(x^{(k)})$ and set $x^{(k+1)} = x^{(k)} + \delta$.

(a) Derive the Jacobian

Write out the Jacobian $J(x) = \partial f_i / \partial x_j$ by hand in the markdown cell below (double-click to edit; use LaTeX). Getting this on paper first is the whole game — part (c) will show you why.

TODO (your derivation here)

$$J(x_1, x_2) = \begin{pmatrix} ? & ? \\ ? & ? \end{pmatrix}$$

(b) Implement Newton's method

Implement the iteration. Requirements:

• Solve $J\delta = -f$ with np.linalg.solve — do not form $J^{-1}$ with np.linalg.inv.
• Start from $x^{(0)} = (0.1,\, 0.1)$.
• Stop when $\lVert f(x^{(k)})\rVert < 10^{-12}$ or after 50 iterations.
• Record $\lVert f(x^{(k)})\rVert$ at every iteration for part (c).

Print the root and the number of iterations.

Cross-check (allowed here): scipy.optimize.root should agree with your root to machine precision.

def f(x):
    return np.array([
        5 * x[0] ** 2 - x[1] ** 2,
        x[1] - 0.25 * (np.sin(x[0]) + np.cos(x[1])),
    ])


def jacobian(x):
    # TODO: fill in the four entries you derived in (a)
    ...


def newton(f, jac, x0, tol=1e-12, max_iter=50):
    # TODO: implement the iteration; return the root and the list of
    # residual norms ||f(x_k)|| (one entry per iteration, including x_0).
    ...


# TODO: run from x0 = (0.1, 0.1); print root and iteration count.
...

(c) Certify it: the verification contract

Your code from (b) prints a root. So would a subtly wrong implementation — this is the central fact of numerical computing (and of working with AI assistants). Two checks certify your work; each catches a classic mistake that produces a perfectly plausible-looking answer.

C1 — Jacobian vs. finite differences. At 5 random points $x \in [-1,1]^2$, compare every entry of jacobian(x) against a central finite difference of f: $\;J_{ij} \approx \bigl(f_i(x + h e_j) - f_i(x - h e_j)\bigr)/2h$ with $h = 10^{-6}$. Report the maximum absolute error over all points and entries — it should be $\sim 10^{-7}$ or smaller. A sign error in one entry shows up as an error the size of the entry itself.

C2 — Convergence rate. Quadratic convergence means $r_{k+1} \approx C\, r_k^2$ where $r_k = \lVert f(x^{(k)})\rVert$: the ratio $r_{k+1}/r_k^2$ stays roughly constant while $r_{k+1}/r_k$ collapses. Plot $r_k$ on a log axis and print both ratio columns. A merely-linear iteration (constant $r_{k+1}/r_k$) converges too — the rate is what certifies your Jacobian and your update.

def fd_jacobian(f, x, h=1e-6):
    n = len(x)
    J = np.zeros((n, n))
    for j in range(n):
        e = np.zeros(n)
        e[j] = h
        J[:, j] = (f(x + e) - f(x - e)) / (2 * h)
    return J


# TODO C1: max |jacobian(x) - fd_jacobian(f, x)| over 5 random points
# (use rng = np.random.default_rng(0)), then assert it is < 1e-6.
...

# TODO C2: semilogy plot of residuals from (b), and a table of
# k, r_k, r_{k+1}/r_k, r_{k+1}/r_k^2.
...

TODO (what does the rate table show, and why is the rate — not convergence — the test?)

(d) A start that cannot work

The origin looks like the most natural starting guess of all. Try $x^{(0)} = (0, 0)$. What happens, and why? Look at $J(0,0)$ before you blame your code.

# TODO: call your newton() from (0, 0) and inspect J(0,0).
...

TODO (your answer here)

Exercise 3 — Newton's method for $\sqrt{34}$

Compute $\sqrt{34}$ with your own scalar Newton iteration: formulate the problem as finding the root of $f(x) = x^2 - 34$, start from $x^{(0)} = 6$, and iterate until $|f(x^{(k)})| < 10^{-12}$. Compare with np.sqrt(34).

# TODO: scalar Newton for f(x) = x^2 - 34, f'(x) = 2x, from x0 = 6.
...