Model fitting

By Anqi Wu
 Who is Anqi Wu
- Postdoc at Columbia with Liam Paninski and John Cunningham

- PhD at Princeton with Jonathan Pillow

- Research interest: neural sensory encoding, latent variable models for

large-scale neural recordings, fMRI decoding, behavior analysis

- Modeling: probabilistic graphical model, Gaussian process, latent

variable/dynamic model, Bayesian deep learning, Bayesian optimization
and active learning

Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 1 2

 Who is Pierre-Etienne
- PhD candidate in Eero Simoncelli’s lab at NYU

- MSc with Sophie Denève and Christian Machens at ENS Paris

- Research interests:

Principles of biological information processing (gasp!)

Visual motion perception and temporal prediction

- Tools: signal processing, statistics and machine learning

Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 1 3

 D3 Team

Alex Hyafil Jan Pierre Etienne Anqi Wu

Kunlin Wei
Drugowitsch Fiquet

Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 1 4

 Roadmap of Week 1 Day 3
Tutorial part 1 Tutorial part 2
- Linear regression - Bias-variance tradeoff
- Mean squared error - Underfitting and overfitting
- Maximum likelihood estimation - Cross validation
- Bootstrapping - K-fold
- Multiple linear regression - Leave-one-out
- Polynomial regression

Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 1 5

 W1D3 Tutorial 1:
Linear Regression &
Mean Squared Error

Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 1 6

 Linear regression makes predictions about the linear relationship between the input
variable (contrast) and the output variable y (neural response).

y = ✓1 x + ✓0
neural response linear weight contrast Intercept

We are not considering the intercept for simplicity, resulting

in a one-parameter model.

Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 1 7

 Linear regression makes predictions about the linear relationship between the input
variable (contrast) and the output variable y (neural response).

y = ✓ x
x
y =✓
neural response linear weight contrast
?
We are not considering the intercept for simplicity, resulting
in a one-parameter model.

Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 1 8

 Mean Squared Error (MSE)
MSE computes the average error between the model prediction ŷ and the true y.

model prediction residual

N
1 X 1 XN
MSE = (yi ŷi )2 = e2i yi
ei {
N i=1 N i=1 ŷi

total number of true neural index of data points

response i=1,…,N
data points

xi

Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 1 9

 Mean Squared Error (MSE)
N
1 X
N
1 X
MSE = (yi 2
ŷi ) = (yi ✓xi )2
N i=1 N i=1
N
N
X
⇤ 1 X
Optimal ✓ =1 argmin
N
X MSE =1 argmin
N
X = (✓x yii))22
(yi i ✓x
✓⇤ = argmin (✓x✓i⇤ =yiargmin
)2 (✓xi yN 2
i ) i=1 i=1
N i=1
N i=1
✓MSE
To minimize MSE, we solve for where its gradient is 0: ✓true
N NX NN N N N
1 X 12 X
N
@MSE 2 2X
12 X
X 2 2
X X
=r✓rr
MSE (y=i =
✓✓MSE
MSE =
=✓x =i )(✓x (y
i i i(y
(✓x
(✓x iyii✓x
)x i✓x
)=i ii=
yyi )x
)x 0)=00 ) ✓ x2i x i yi = 0
N@✓ NN NN
i=1
N
i=1 i=1
i=1i=1
i=1 i=1 i=1
PN
x i yi
) ✓MSE = Pi=1
N
x2i ✓
i=1

Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 1 10

 Notebook
In the following notebook, you will
1. Define a linear model and generate some simulated data from this underlying linear model.
2. Fit linear regression to the simulated data by minimizing MSE, and recover the linear weight.

Enjoy!!

Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 1 11

 W1D3 Tutorial 2:
Maximum Likelihood
Estimation

Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 2 12

 Maximum Likelihood Estimation (MLE)
Gaussian/Normal distribution
Data can be "distributed" (spread out) in different ways. But there are many cases where the data tends to be
around a central value with no bias left or right, and it gets close to a "Gaussian Distribution" like this:

Bell Curve
A Gaussian distribution function is defined by

mean parameter µ standard deviation >0

2
variance f (x)
= 15
1 1
(x µ)2
probability density function f (x) = p e 2 2
µ = 100
2⇡2⇡
x
Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 2 13
 Maximum Likelihood Estimation (MLE)

If we generate y from x using However this is the data we observe, with noise.
y = ✓ x

all data on the line what we don’t care about

without noise (e.g., deviation from mean firing rate)
what we can’t control
(e.g., measurement noise)

Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 2 14

 Maximum Likelihood Estimation (MLE)

If we generate y from x using Model data with noise.

y = ✓ x y = ✓ x + ⌘ Noise
all data on the line what we don’t care about
without noise (e.g., deviation from mean firing rate)
what we can’t control
(e.g., measurement noise)
2
⌘ ⇠ N (0, )

2
y ⇠ N (✓x, )

Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 2 15

 Maximum Likelihood Estimation (MLE)
We consider the linear regression model with a Gaussian distributed noise ⌘ .

y = ✓ x + ⌘ encoding model
1 1
(y ✓x)2
p(y|x, ✓) = p e 2 2

neural response linear weight contrast Gaussian 2⇡2⇡

noise mean(y) = ✓x
2
⌘ ⇠ N (0, ) 2
2
var(y) =
the
the larger
larger is, the noisier the data is
the noisier the data is

Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 2 16

 Maximum Likelihood Estimation (MLE)
1 1
(y ✓x)2
p(y|x, ✓) = e 2 2
2⇡
p(y|x = 8, ✓)
encoding model

1 1
(y 8✓)2
p(y|x = 8, ✓) =p e 2 2
2⇡
2
mean(y) = 8✓
8✓ var(y) = 2

Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 2 17

 Maximum Likelihood Estimation (MLE)
1 1
(y ✓x)2
p(y|x, ✓) = e 2 2
2⇡
p(y|x = 2, ✓)
encoding model

1 1
(y 2✓)2
p(y|x = 2, ✓) =p e 2 2
2⇡
2
mean(y) = 2✓
2✓ 2
var(y) =

Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 2 18

 Maximum Likelihood Estimation (MLE)
Say it here
observed data {x, y}

encoding model
✓✓✓ =
= 0.5
= 2.2
1.2 <latexit sha1_base64="ma+UDoqMeaR0Nol6xlhWbd5k0II=">AAACB3icbVDLSgNBEJyNrxhfUY9eFoPgadmNil6EoBePEcwDkyXMTnqTIbOzy0yvEEI+wLtX/QVv4tXP8A/8DCfJHkxiQUNR1U13V5AIrtF1v63cyura+kZ+s7C1vbO7V9w/qOs4VQxqLBaxagZUg+ASashRQDNRQKNAQCMY3E78xhMozWP5gMME/Ij2JA85o2ikxzb2Aem161x0iiXXcaewl4mXkRLJUO0Uf9rdmKURSGSCat3y3AT9EVXImYBxoZ1qSCgb0B60DJU0Au2PpheP7ROjdO0wVqYk2lP178SIRloPo8B0RhT7etGbiP95rRTDK3/EZZIiSDZbFKbCxtievG93uQKGYmgIZYqbW23Wp4oyNCHNbdGMClAgxgUTjbcYxDKplx3vzCnfn5cqN1lIeXJEjskp8cglqZA7UiU1wogkL+SVvFnP1rv1YX3OWnNWNnNI5mB9/QJJ0Zlq</latexit>

1 (y ✓x)2
<latexit sha1_base64="1HLYojB8bxtIX8VdSbQk1zvcdq0=">AAACB3icbVDLSsNAFJ34rPVVdekmWARXIamCboSiG5cV7APbUCbTm3boZBJmboQS+gHu3eovuBO3foZ/4Gc4bbOwrQcuHM65l3vvCRLBNbrut7Wyura+sVnYKm7v7O7tlw4OGzpOFYM6i0WsWgHVILiEOnIU0EoU0CgQ0AyGtxO/+QRK81g+4CgBP6J9yUPOKBrpsYMDQHpdcSrdUtl13CnsZeLlpExy1Lqln04vZmkEEpmgWrc9N0E/owo5EzAudlINCWVD2oe2oZJGoP1sevHYPjVKzw5jZUqiPVX/TmQ00noUBaYzojjQi95E/M9rpxhe+RmXSYog2WxRmAobY3vyvt3jChiKkSGUKW5utdmAKsrQhDS3RTMqQIEYF0003mIQy6RRcbxzp3J/Ua7e5CEVyDE5IWfEI5ekSu5IjdQJI5K8kFfyZj1b79aH9TlrXbHymSMyB+vrF0g7mWk=</latexit>

<latexit sha1_base64="MNH2uR2LpvXiOeYyDoIltkSnpPA=">AAACB3icbVDLSsNAFJ34rPVVdekmWARXIamCboSiG5cV7APbUCbTm3boZBJmboQS+gHu3eovuBO3foZ/4Gc4bbOwrQcuHM65l3vvCRLBNbrut7Wyura+sVnYKm7v7O7tlw4OGzpOFYM6i0WsWgHVILiEOnIU0EoU0CgQ0AyGtxO/+QRK81g+4CgBP6J9yUPOKBrpsYMDQHrtOZVuqew67hT2MvFyUiY5at3ST6cXszQCiUxQrduem6CfUYWcCRgXO6mGhLIh7UPbUEkj0H42vXhsnxqlZ4exMiXRnqp/JzIaaT2KAtMZURzoRW8i/ue1Uwyv/IzLJEWQbLYoTIWNsT153+5xBQzFyBDKFDe32mxAFWVoQprbohkVoECMiyYabzGIZdKoON65U7m/KFdv8pAK5JickDPikUtSJXekRuqEEUleyCt5s56td+vD+py1rlj5zBGZg/X1C0afmWg=</latexit>

1
p(y|x, ✓) = p e 2 2
2⇡2⇡
mean(y) = ✓x 1 1
(
p(y|x, ✓) is=unknown
e 2 2
var(y) = 2 2⇡

Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 2 19

 Maximum Likelihood Estimation (MLE)
observed data {x, y}

L(✓|x, y) = p(y|x,
Likelihood: which ✓)leads to the Gaussian distribution
that most likely matches the data points {x, y}?

1 1
(y ✓x)2
L(✓|x, y) = p(y|x, ✓)
✓) = p e 2 2
2⇡
2⇡

<latexit sha1_base64="8Om5nnTmDn6YlI40fdeYnVx+XGU=">AAACEHicbVDLSgNBEJyNrxhfGz16WQxCBAm7UdGLEPTiMYJ5QBLC7KQ3GTL7YKZXDTE/4d2r/oI38eof+Ad+hpNkDyaxoKGo6qaaciPBFdr2t5FaWl5ZXUuvZzY2t7Z3zOxuVYWxZFBhoQhl3aUKBA+gghwF1CMJ1HcF1Nz+9div3YNUPAzucBBBy6fdgHucUdRS28xG+cHT43ETe4D00i6cHbXNnF2wJ7AWiZOQHElQbps/zU7IYh8CZIIq1XDsCFtDKpEzAaNMM1YQUdanXWhoGlAfVGs4eX1kHWqlY3mh1BOgNVH/Xgypr9TAd/WmT7Gn5r2x+J/XiNG7aA15EMUIAZsGebGwMLTGPVgdLoGhGGhCmeT6V4v1qKQMdVszKYpRARLEKKOrceaLWCTVYsE5KRRvT3Olq6SkNNknByRPHHJOSuSGlEmFMPJAXsgreTOejXfjw/icrqaM5GaPzMD4+gWtY5w7</latexit>
p(y|x, ✓ = 0.5) <latexit sha1_base64="hrKwlJcPmCMShMS67F2YIkkBldc=">AAACEHicbVDLSgNBEJz1GeMr0aOXwSBEkLAbBb0IQS8eI5gHJCHMTjrJkNkHM73qsuYnvHvVX/AmXv0D/8DPcPI4mMSChqKqm2rKDaXQaNvf1tLyyuraemojvbm1vbObye5VdRApDhUeyEDVXaZBCh8qKFBCPVTAPFdCzR1cj/zaPSgtAv8O4xBaHuv5ois4QyO1M9kwHz89njSxD8gunULxuJ3J2QV7DLpInCnJkSnK7cxPsxPwyAMfuWRaNxw7xFbCFAouYZhuRhpCxgesBw1DfeaBbiXj14f0yCgd2g2UGR/pWP17kTBP69hzzabHsK/nvZH4n9eIsHvRSoQfRgg+nwR1I0kxoKMeaEco4ChjQxhXwvxKeZ8pxtG0NZOiOZOgQA7TphpnvohFUi0WnNNC8fYsV7qalpQiB+SQ5IlDzkmJ3JAyqRBOHsgLeSVv1rP1bn1Yn5PVJWt6s09mYH39AqovnDk=</latexit>
p(y|x, ✓ = 1.2) <latexit sha1_base64="vEb7HfjhUP/Ma6NQU4iwRwp8m3A=">AAACEHicbVDLSgNBEJz1GeNro0cvi0GIIMvuKuhFCHrxGME8IAlhdtJJhsw+mOlVl5if8O5Vf8GbePUP/AM/w8njYBILGoqqbqopPxZcoeN8G0vLK6tr65mN7ObW9s6umdurqCiRDMosEpGs+VSB4CGUkaOAWiyBBr6Aqt+/HvnVe5CKR+EdpjE0A9oNeYczilpqmbm4kD49njSwB0gvPds7bpl5x3bGsBaJOyV5MkWpZf402hFLAgiRCapU3XVibA6oRM4EDLONREFMWZ92oa5pSANQzcH49aF1pJW21YmknhCtsfr3YkADpdLA15sBxZ6a90bif149wc5Fc8DDOEEI2SSokwgLI2vUg9XmEhiKVBPKJNe/WqxHJWWo25pJUYwKkCCGWV2NO1/EIql4tntqe7dn+eLVtKQMOSCHpEBcck6K5IaUSJkw8kBeyCt5M56Nd+PD+JysLhnTm30yA+PrF6vMnDo=</latexit>
p(y|x, ✓ = 2.2)
✓ = 0.5
<latexit sha1_base64="ma+UDoqMeaR0Nol6xlhWbd5k0II=">AAACB3icbVDLSgNBEJyNrxhfUY9eFoPgadmNil6EoBePEcwDkyXMTnqTIbOzy0yvEEI+wLtX/QVv4tXP8A/8DCfJHkxiQUNR1U13V5AIrtF1v63cyura+kZ+s7C1vbO7V9w/qOs4VQxqLBaxagZUg+ASashRQDNRQKNAQCMY3E78xhMozWP5gMME/Ij2JA85o2ikxzb2Aem161x0iiXXcaewl4mXkRLJUO0Uf9rdmKURSGSCat3y3AT9EVXImYBxoZ1qSCgb0B60DJU0Au2PpheP7ROjdO0wVqYk2lP178SIRloPo8B0RhT7etGbiP95rRTDK3/EZZIiSDZbFKbCxtievG93uQKGYmgIZYqbW23Wp4oyNCHNbdGMClAgxgUTjbcYxDKplx3vzCnfn5cqN1lIeXJEjskp8cglqZA7UiU1wogkL+SVvFnP1rv1YX3OWnNWNnNI5mB9/QJJ0Zlq</latexit>

✓ = 1.2
<latexit sha1_base64="MNH2uR2LpvXiOeYyDoIltkSnpPA=">AAACB3icbVDLSsNAFJ34rPVVdekmWARXIamCboSiG5cV7APbUCbTm3boZBJmboQS+gHu3eovuBO3foZ/4Gc4bbOwrQcuHM65l3vvCRLBNbrut7Wyura+sVnYKm7v7O7tlw4OGzpOFYM6i0WsWgHVILiEOnIU0EoU0CgQ0AyGtxO/+QRK81g+4CgBP6J9yUPOKBrpsYMDQHrtOZVuqew67hT2MvFyUiY5at3ST6cXszQCiUxQrduem6CfUYWcCRgXO6mGhLIh7UPbUEkj0H42vXhsnxqlZ4exMiXRnqp/JzIaaT2KAtMZURzoRW8i/ue1Uwyv/IzLJEWQbLYoTIWNsT153+5xBQzFyBDKFDe32mxAFWVoQprbohkVoECMiyYabzGIZdKoON65U7m/KFdv8pAK5JickDPikUtSJXekRuqEEUleyCt5s56td+vD+py1rlj5zBGZg/X1C0afmWg=</latexit>
✓ = 2.2
<latexit sha1_base64="1HLYojB8bxtIX8VdSbQk1zvcdq0=">AAACB3icbVDLSsNAFJ34rPVVdekmWARXIamCboSiG5cV7APbUCbTm3boZBJmboQS+gHu3eovuBO3foZ/4Gc4bbOwrQcuHM65l3vvCRLBNbrut7Wyura+sVnYKm7v7O7tlw4OGzpOFYM6i0WsWgHVILiEOnIU0EoU0CgQ0AyGtxO/+QRK81g+4CgBP6J9yUPOKBrpsYMDQHpdcSrdUtl13CnsZeLlpExy1Lqln04vZmkEEpmgWrc9N0E/owo5EzAudlINCWVD2oe2oZJGoP1sevHYPjVKzw5jZUqiPVX/TmQ00noUBaYzojjQi95E/M9rpxhe+RmXSYog2WxRmAobY3vyvt3jChiKkSGUKW5utdmAKsrQhDS3RTMqQIEYF0003mIQy6RRcbxzp3J/Ua7e5CEVyDE5IWfEI5ekSu5IjdQJI5K8kFfyZj1b79aH9TlrXbHymSMyB+vrF0g7mWk=</latexit>
L(✓ = 1.2) > L(✓ = 0.5) > L(✓ = 2.2)
<latexit sha1_base64="3Hq/vSdbFpcbXXwyfbfAc6MfcSg=">AAACRnicbVC7SgNBFL0b3/EVtbRZDEJslt2oaKMELbSwUDAPSEKYndwkQ2YfzNwVwpJP8j/srQTtrezE1klMoSYHBg7nnMu9c/xYCk2u+2Jl5uYXFpeWV7Kra+sbm7mt7YqOEsWxzCMZqZrPNEoRYpkESazFClngS6z6/cuRX31ApUUU3tMgxmbAuqHoCM7ISK3cVSNg1ONMpjfDQoN6SOzMc4oH5zN01zmeqRdNvpXLu447hj1NvAnJwwS3rdx7ox3xJMCQuGRa1z03pmbKFAkucZhtJBpjxvusi3VDQxagbqbjDw/tfaO07U6kzAvJHqu/J1IWaD0IfJMcnav/eyNxlldPqHPaTEUYJ4Qh/1nUSaRNkT1qz24LhZzkwBDGlTC32rzHFONkOv6zRZuOUKEcZk013v8ipkml6HiHTvHuKF+6mJS0DLuwBwXw4ARKcA23UAYOj/AMr/BmPVkf1qf19RPNWJOZHfiDDHwDUCCwXA==</latexit>

L(✓|x, y)
Therefore, we could calculate the=optimal
p(y|x, ✓)
by maximizing the likelihood.

Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 2 20

 Maximum Likelihood Estimation (MLE)
Empirically, we maximize the log likelihood
N
Y XN N
1 11 X
11 1
(yi ✓xi )2 22
log L(✓|x,y)
log L(✓|x, y) = log
= log
log p 22 e 2 2(y
(y i ✓xii))
✓x
2⇡ 222⇡N
2⇡ N i=1 i
i=1 i=1
N
X 1 1
(yi ✓xi )2
= log p e 2 2

i=1
2⇡
N
X X N
1 1
(yi ✓xi )2
= log p + log e 2 2

i=1
2⇡ i=1
X NX
N
N1 1 21 1 1 2 2
L(✓|x,
loglog L(✓|x,
y)y)
==logloglog 2⇡ 2 2 2 (yi(yi ✓x✓x
i) i)
2 2⇡2⇡ 2 2 N2N
i=1
i=1

Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 2 21

 Maximum Likelihood Estimation (MLE)
N N
1 X1 X
Empirically, we maximize the log likelihood = = (ŷi (ŷyi )2 yi )2
Nx MSEMSE
N i=1
N i=1
N
NXN
N 11 1 11 1 X
2 1X 22 2
log L(✓|x,
L(✓|x,
loglog L(✓|x,
y)
y)y)
== log
=logloglog 2⇡ 22 2 2 (y
(yii(yi ✓x
✓xi✓x
i)) i ) ✓MLE
22⇡
2⇡2⇡ 22 2NN2N i=1
i=1
i=1

Differentiate and set to zero ✓MSE

N ✓true
@log L(✓|x, y) 11 X
r✓ log L(✓|x, y) == 2 (yi ✓xi )xi = 0
@✓ N2 i=1
PN
x i yi
) ✓MLE = Pi=1
N
i=1 x2i
PN
i=1 xi yi ✓
) ✓MSE =
(the same as when
PNp(y|x, ✓) is Gaussian)
x 2
i=1 i
Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 2 22
 Notebook
In the following notebook, you will
1. Code for likelihood function for Gaussian distribution.
2. Fit linear regression with maximum likelihood estimation.

Enjoy!!

Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 2 23

 W1D3 Tutorial 3:
Bootstrapping

Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 3 24

 Uncertainty of Parameters

MSE = 0.38

Question: How confident we are about the model fit?

Use resampling to estimate uncertainty thus confidence

Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 3 25

 Resampling to Estimate Uncertainty
Knowing the true distribution

0.9 1.0 1.1 1.2

Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 3 26

 Resampling to Estimate Uncertainty
Knowing the underlying true distribution
95% of bootstrap
samples ➝
confidence intervals

✓est

…
0.9 1.0 1.1 1.2

Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 3 27

 Assessing Uncertainty by Bootstrap
Of course we never know the true population; we only have data samples.

data distribution observed data samples

Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 3 28

 Assessing Uncertainty by Bootstrap
1. Resampling from the observed dataset with replacement; the number of the
resampled data points is the same as the number of observed data points.
2. Collect all estimates into a distribution, and analyze the confidence intervals.
95% of bootstrap
samples ➝
confidence intervals
✓true ✓est

…
0.9 1.0 1.1 1.2

bad estimate
Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 3 29
 Assessing Uncertainty by Bootstrap
1. Resampling from the observed dataset.
2. Collect all estimates into a distribution, and analyze the confidence intervals.
95% of bootstrap
samples ➝
confidence intervals
✓true ✓est

0.9 1.0 1.1 1.2

good estimate
Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 3 30
 Assessing Uncertainty by Bootstrap

95% of bootstrap
samples ➝
confidence intervals • In most of real-world applications, we don’t know the
✓true ✓est ground truth.
• But it’s still beneficial to have a distribution rather
than a point estimate.
• With distribution and uncertainty, you have more
0.9 1.0 1.1 1.2
information to make decision.

Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 3 31

 Notebook
In the following notebook, you will
1. Generate simulated dataset from an underlying linear model.
2. Implement bootstrap to get the confidence interval of the model fit.
3. Compare the true parameter with the estimated distribution to see how confident the model is.

Enjoy!!

Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 3 32

 W1D3 Tutorial 4:
Multiple
Linear Regression

Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 4 33

 Linear model y = ✓ x

neural response linear weight contrast

Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 4 34

 [Simple] y y= ✓=
0 + ✓1 x1yx+=✓2✓x
0y2++
=✓...
1✓x
01++✓✓d✓1x2xdx1 2+
+ ++⌘✓...
2 x+ ✓d
2+
Linear model
neural response linear weight contrast intercept

Multiple y = ✓0 + ✓1 x1 + ✓2 x2 + ... + ✓d xd + ⌘
linear model
neural response Assume multiple stimulus features
(e.g., orientation, contrast, etc.)

Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 4 35

 More generally, >
in vector form,
y i = ✓ xi 8i = 1, ..., N

for one data point i >

✓ = [✓0 , ✓1 , ✓2 ,x...]=>[1, xi,1
=, x ...], >
[1,i,2x,i,1 xi,2 , ...]

Build up to y = X✓
matrix version
y1 y2 y0
Index i x1 x2 x0 ✓ 1 ✓ 2 ✓ 0
y1 y2 y0 x1 x2 x0 ✓ 1 ✓ 2 ✓ 0
y1 y2 y0y = Xk + ⌘
x1 x2 x0 ✓ 1 ✓ 2 ✓ 0

design matrix
Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 4 36
 Geometric Interpretation of MSE
MSE computes the average error between the model prediction ŷ and the true y.

model prediction residual

N
1 X 1 XN
MSE = (yi ŷi )2 = e2i yi
ei {
N i=1 N i=1 ŷi

total number of true neural index of data points

response i=1,…,N
data points

xi

Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 4 37

 Geometric Interpretation of MSE

2D case without intercept e=y (✓1 x1 + ✓2yx2=) ✓ x + ✓ x

1 1 2 2

model predict ŷ = ✓1 x1 + ✓2 x2 y = ✓1 x1
Residual + ✓2 x2
y = ✓1 x1 + ✓2 x2
residual e=y ŷ(✓
e= 1 x1
y + ✓(✓2 x
1x2 )1 + ✓2 x2 ) ŷ = ✓1 x1 + ✓2 x2
y = ✓1 x1 + ✓2 x2 Regression Plane

minimize MSE = (✓1 x1 + ✓2 x2 y)2 = e2

y = ✓1 x1 + ✓2 x2

source: [Link]

Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 4 38

 the number of data points
MSE solution:
XN
⇥ 2 > ⇤
⇤
✓ = argmin ||X✓
Tr (X✓ y)>> (X✓
⇥ y||2 = (yi ⇤ ✓y)xi )2
⇤
✓ = argmin Tr (X✓ y) i=1
(X✓ y)

Differentiate and set to zero y = X✓

⇥ >
⇤
@ Tr (X✓ @y)||X✓ y)22
(X✓ y||
= 2X > (X✓ y) = 0
@✓ @✓
y1 y2 y0 x1 x2 x0 ✓ 1 ✓ 2 ✓ 0

Index i
y1 y2 y0 x1 x2 x0 ✓ 1 ✓ 2 ✓ 0
(X>>X)
X) 1 1XX>>yy
) ✓
kMSE = (X
MSE =
y = Xk + ⌘
y1 y2 y0 x1 x2 x0 ✓ 1 ✓ 2 ✓ 0

design matrix
Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 4 39
 MLE solution:
YN N
X
d 2 1 >
argmax log L(✓|X, y) = log log 2⇡ ✓) = log
p(y|X, (y X✓) |xi , ✓)
p(yi(y X✓)= log p(yi |xi , ✓)
2 2 2
i=1 i=1

=
N
X
log p
1
e
1
2 2
(yi ✓ > xi )2 yy =
= X✓
Xk + ⌘
i=1
( 2⇡)d the number of
input features
N
X X N
1 1
(yi ✓ > xi )2 y1 y2 y0 x1 x2 x0 ✓ 1 ✓ 2 ✓ 0 ⌘ 1 ⌘ 2 ⌘ 0

Index i
= log p + log e 2 2

✓1 ✓2 ✓0 ⌘ 1 ⌘ 2 ⌘ 0
( 2⇡) d
i=1 i=1 y1 y2 y0 x1 x2 x
y 0= Xk +⌘
Nd
N
1 X
y = Xk + ⌘
y1 y2 y0 x1 x2 x0 ✓ 1 ✓ 2 ✓ 0 ⌘ 1 ⌘ 2 ⌘ 0
= log 2⇡ 2
2
(yi ✓ > xi ) 2
2 2 i=1
XN 2
N dd 11 > >
⌘ ⇠ N (0, )
222 2
g L(✓|X,=
y) = log log 2⇡
2⇡
2⇡ 22
(y (y
X✓)
i ✓(y xi )X✓)
22 2
i=1 design matrix
Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 4 40
 MLE solution:
N
⇤ d 2 1 N dd > 222 11 X > > 2
argmax log L(✓|X, = log log
✓ =y)argmax p(y|X, ✓) =
log2⇡L(✓|X, (y log
y)2 = log 2⇡
X✓) 2⇡(y X✓)22(y (y
2⇡ X✓)
i ✓(y xi )X✓)
2 2 22 2
i=1

Differentiate and set to zero

yy =
= X✓
Xk + ⌘
@ log L(✓|X, y) 1
= X > (X✓ y) = 0
y1 y2 y0 x1 x2 x0 ✓ 1 ✓ 2 ✓ 0 ⌘ 1 ⌘ 2 ⌘ 0

Index i
@✓ 2

y1 y2 y0 x1 x2 x ✓1 ✓2 ✓0 ⌘ 1 ⌘ 2 ⌘ 0
(X>>X)
X) 1 1XX>>yy y 0= Xk +⌘
) ✓✓
kMLE = (X
MSE =
MSE
y = Xk + ⌘
y1 y2 y0 x1 x2 x0 ✓ 1 ✓ 2 ✓ 0 ⌘ 1 ⌘ 2 ⌘ 0
(same as MSE when the noise is Gaussian)
2
⌘ ⇠ N (0, )

design matrix
Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 4 41
 Polynomial Regression
P
X
y = ✓0 + ✓1 x + ✓2 x2 + ... + ✓P xP = ✓p xp
p=0
More generally, Build up to
in vector form, matrix version
y = X✓
for one data point i y 1 y 2 y 0 x1 x2 x0

Index i
✓1 ✓2 ✓0
y1 y2 y0 x1 x2 x0 ✓ 1 ✓ 2 ✓ 0
y i = ✓ > xi 8i = 1, ..., N y1 y2 y0 y = Xk + ⌘
x1 x2 x0 ✓ 1 ✓ 2 ✓ 0
The only difference

✓ = [✓0 , ✓1 , ✓2 , ...]> xi = [1, xi , x2i , ...]>

=, x
x = [1, xi,1 [1,i,2x,i,1 xi,2 , ...]>
...], > linear design matrix
Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 4 42
 Notebook
In the following notebook, you will
1. Build design matrix for multiple linear regression.
2. Calculate MSE/MLE solutions with matrix forms.
Optional exercise:
1. Generate simulated dataset from a polynomial model.
2. Build design matrix with increasing powers. Enjoy!!
3. Calculate MSE/MLE solutions with matrix forms.

Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 4 43

 W1D3 Tutorial 5:
Bias-variance Tradeoff

Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 5 44

First order polynomial Higher order
(Linear regression) Second order polynomial polynomial
5
X
y = ✓1 x + ✓0 y = ✓2 x2 + ✓1 x + ✓0 y= ✓p xp
<latexit sha1_base64="M3lPwR+tLR+20fGNiDiqWPB9Oig=">AAACGHicbVDJSgNBEO1xjXGLehIvg0HwFGaiopdA0IvHCGaBLENPp5I06ZlpumvEMAT/w7tX/QVv4tWbf+Bn2FkOJvFBweO9Kqrq+VJwjY7zbS0tr6yurac20ptb2zu7mb39io5ixaDMIhGpmk81CB5CGTkKqEkFNPAFVP3+zcivPoDSPArvcSChGdBuyDucUTSSlzkcFBo6DrxEFpxh66KBPUDqyceW9DJZJ+eMYS8Sd0qyZIqSl/lptCMWBxAiE1TruutIbCZUIWcChulGrEFS1qddqBsa0gB0Mxm/MLRPjNK2O5EyFaI9Vv9OJDTQehD4pjOg2NPz3kj8z6vH2LlqJjyUMULIJos6sbAxskd52G2ugKEYGEKZ4uZWm/WoogxNajNbNKMCFIhh2kTjzgexSCr5nHuWy9+dZ4vX05BS5Igck1PikktSJLekRMqEkSfyQl7Jm/VsvVsf1uekdcmazhyQGVhfv6idoKQ=</latexit>
p=0

2nd order
polynomial
Which one is
linear the best?

5th order
polynomial

Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 5 45

 Bias and Variance
Bias: systematic deviation from structure underlying data, low model complexity.
Variance: capturing variability in predictions from models trained on different datasets, high model complexity.

True model:
2nd order Low model complexity

bias { variance Bias is large: the linear slope always

misses many data points
linear Variance is small: linear slopes are
close to each other across samples

Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 5 46

 Bias and Variance
Bias: systematic deviation from structure underlying data, low model complexity.
Variance: capturing variability in predictions from models trained on different datasets, high model complexity.

{
5th order True model:
polynomial 2nd order High model complexity
Bias is small: the 5th order polynomials go
through almost all data points
variance
bias Variance is large: the 5th order polynomials
are much different from each other

Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 5 47

 Bias-variance Trade-off
Total error = bias + variance Best model: balances bias / variance

best
2nd order model
polynomial

model error
linear

5th order total

bias variance
polynomial

model complexity
underfitting overfitting
Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 5 48
 Bias-variance Trade-off
Model evaluation aims at estimating the generalization error of the selected model, i.e.,
how well the selected model performs on unseen data.
Training set: data used for fitting Test set: unseen data (not used for fitting)
best test set
2nd order model error (total)
polynomial

linear
test set training set
error (bias)
5th order
polynomial
training set

Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 5 49

 Notebook
In the following notebook, you will
1. Continue the polynomial regression exercise.
2. Split the data into training and test sets.
3. Fit different polynomial models using the training set and test on the test set.
4. Inspect the training error and the test error.
5. Explore overfitting and underfitting.

Enjoy!!

Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 5 50

 W1D3 Tutorial 6:
Cross Validation
We can select the right complexity model in a
data-driven way.

Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 6 51

 K-fold
Cross Cross
Validation
Validation
Motivation: one random split could result in a biased training-test distribution
Solution: multiple random splits with k folds

Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 6 52

 K-fold Cross Validation
Motivation: one random split could result in a biased training-test distribution
Solution: multiple random splits with k folds.

Steps: P
X
0. Define a model to evaluate, e.g. y= ✓p xp when P=1
p=0

1. Take one fold as the test set, the rest as the training set.
2. Fit the model with the training set, and get the optimal ✓p .
3. Estimate MSE on the test set using the above ✓p ,
N
X test P
X
1 2
MSEtest = ( ✓p xptest ytest )
Ntest i=1 p=0

where {xtest , ytest }N

i=1 denotes the test set.
test

Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 6 53

 K-fold Cross Validation
Motivation: one random split could result in a biased training-test distribution
Solution: multiple random splits with k folds.

model
P
X
y= ✓p xp
p=0

when P=1

MSE-fold1 MSE-fold2 MSE-fold3 MSE-fold4

averaged MSE
Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 6 54
 K-fold Cross Validation
Motivation: one random split could result in a biased training-test distribution
Solution: multiple random splits with k folds.
best
model

model
P
X Repeat the same 4-fold cross validation
y= ✓p xp
on other Ps, i.e. P = {2, 3, 4…}, get
p=0
averaged test MSE
when P=1
averaged MSE

Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 6 55

 Leave-one-out Cross Validation
Instead of leaving a fold out for test, leaving one data point out.

k-fold leave-one-out

Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 6 56

 Leave-one-out Cross Validation
Instead of leaving a fold out for test, leaving one data point out.

…
best
model

N times

Repeat the same leave-one-

averaged MSE out cross validation on all Ps,
i.e. P = {1, 2, 3, 4…}, get
averaged test MSE

Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 6 57

 Notebook
In the following notebook, you will
1. Implement K-fold cross validation for model selection of polynomial regression.

Enjoy!!

Anqi Wu ⦁ Model Fitting Week 1 ⦁ Day 3 ⦁ Tutorial 6 58