M7777 Applied Functional Data Analysis 4. From Data to Functions — Smoothing Penalties Jan Koláček (kolacek@math.muni.cz) Dept. of Mathematics and Statistics, Faculty of Science, Masaryk University, Brno n Koláček (SCI MUNI) M7777 Applied FDA Fall 2019 1/27 Smoothing Penalties Some disadvantages of basis expansions • Discrete choice of number of basis functions =4> Large effect on results. • Non-hierarchical bases (e.g. B-splines) make life more complicated. Alternatives • Kernel methods A kernel function K[t) at each data point gives weights to observations Use the bandwidth h to regulate smoothness. Jan Koláček (SCI MUNI) M7777 Applied FDA Fall 2019 2 / 27 Smoothing Penalties What do we mean by smoothness? • Some things are fairly clearly smooth: • constants • straight lines • What we really want to do is eliminate small "wiggles" in the data while retaining the right shape. Too smooth Too rough Just right st- Jonns St. Johns St. Johns Jan Koláček (SCI MUNI) M7777 Applied FDA Fall 2019 3 / 27 Smoothing Penalties The D Operator We use the notation that for a function x(t) Dx(t) = -x(t) We can also define further derivatives in terms of powers of D d2 dk °2*(0 = -ro*{t),..., Dkx{t) = -rrx{t),... dt2 dtk • Dx(t) ... the slope of x(t) • D2x(t) ... the curvature of x(t) We measure the size of the curvature for all of x by J2(x) = J [D2x(t)]2dt. Jan Koláček (SCI MUNI) M7777 Applied FDA Fall 2019 4 / 27 Smoothing Penalties Curvature for 3 Fourier Bases Jan Koláček (SCI MUNI) M7777 Applied FDA Fall 2019 5 / 27 Smoothing Penalties Curvature for Fourier Bases Jan Koláček (SCI MUNI) M7777 Applied FDA Fall 2019 6 / 27 Smoothing Penalties Curvature for 35 Fourier Bases Jan Koláček (SCI MUNI) M7777 Applied FDA Fall 2019 7 / 27 Smoothing Penalties Curvature for 75 Fourier Bases Jan Koláček (SCI MUNI) M7777 Applied FDA Smoothing Penalties Penalized Squared Error • We now have two competing desires: fit to data and smoothness. • We will explicitly trade them off by minimizing the penalized sum of squared errors N PENSSEx(x) = ^{yi ~ x(t/))2 + AJ2(x) i=i • A is a smoothing parameter ... measures compromise between fit and smoothness • A —)► oo: roughness increasingly penalized, x(t) becomes smooth. • A —)► 0: penalty reduces, x(t) fits data better. Jan Koláček (SCI MUNI) M7777 Applied FDA Fall 2019 9 / 27 Smoothing Penalties Smoothing with In A St. Johns Days Jan Koláček (SCI MUNI) M7777 Applied FDA Fall 2019 10 / 27 Smoothing Penalties Smoothing with In A St. Johns • • t • • • • • • • • »/L « • • • • • / T * ft "wW • • • • •A* t7tJ • • A •• /\ A • • - • • • * • • • • 0 100 200 300 Days Jan Koláček (SCI MUNI) M7777 Applied FDA Fall 2019 11 / 27 Smoothing Penalties Smoothing with In A St. Johns • • t • • • • • • • # • m • • • » — • • • •• • 4 • • • If Vhf», • • • • • * • Al. • • • • 0 100 200 300 Days Jan Koláček (SCI MUNI) M7777 Applied FDA Fall 2019 12 / 27 Smoothing Penalties Smoothing with In A St. Johns • • t • • • • • • • # • m • • • • - • • • • • —rg--*=■■ • • • • • • • • • • - • • • 0 100 200 300 Days Jan Koláček (SCI MUNI) M7777 Applied FDA Fall 2019 13 / 27 Smoothing Penalties Smoothing with In A St. Johns • • t • • • • • • • # • ( 1 • • • • • •• • • WV v.\: • ___-- • • * On? to • • - • • • • < • • • » ^ « • 0 100 200 300 Days Jan Koláček (SCI MUNI) M7777 Applied FDA Fall 2019 14 / 27 Smoothing Penalties Smoothing with In A St. Johns • • • • • • • • • • • • • • • • •• • 4 • • • •••• • A. ~ w + • • J • • • • • • • • • • • * ■ft * • • 0 100 200 300 Days Jan Koláček (SCI MUNI) M7777 Applied FDA Fall 2019 15 / 27 Smoothing Penalties The Smoothing Spline Theorem • A remarkable theorem tells us that the function x(t) that minimizes N PENSSEx(x) = ^{yi ~ x(ti))2 + AJ2(x) i=i is • a spline function of order 4 (piecewise cubic) • with a knot at each sample point tj • This is often referred to simply as cubic spline smoothing. • The theorem tells us that x(t) takes the form x(t) = 0*(t)c. where 0*(t) is a vector of B-spline basis functions. • The number of basis functions is (A/ — 2) + 4 = N + 2 where N is the number of sampling points. • How do we calculate c? Jan Koláček (SCI MUNI) M7777 Applied FDA Fall 2019 16 / 27 Smoothing Penalties Calculating the Penalized Fit For x(t) = 0*(t)c we have Jm(x) = J [Dmx(t)]2dt = J c/Dm0*(t)/Dm0*(t)cc/t = c/Rc. R is known as the penalty matrix. The PENSSEx takes the form PENSSEx(x) = (y - 0c)'(y - 0c) + Ac;Rc. It is minimized by c= (0;0 + AR)_1 0;y. This is still a linear smoother: y = 0(0^0+ AR)"10/y = SAy. Jan Koláček (SCI MUNI) M7777 Applied FDA Fall 2019 17 / 27 Smoothing Penalties Linear Smooths and Degrees of Freedom • In least squares fitting, the degrees of freedom used to smooth the data is exactly K, the number of basis functions. • In penalized smoothing, we can have K > N. • The smoothing penalty reduces the flexibility of the smooth (i.e., we say we know something). • The degrees of freedom are controlled by A. A natural measure turns out to be df(X) = trSA. Jan Koláček (SCI MUNI) M7777 Applied FDA Fall 2019 18 / 27 Smoothing Penalties Degrees of Freedom for Precipitation in St. Johnes Degrees of Freedom 150 H 3 7 l'l 15 ~\9 log lambda Jan Koláček (SCI MUNI) M7777 Applied FDA Fall 2019 19 / 27 Smoothing Penalties Alternative Definitions of Roughness • D2x(t) is only one way to measure the roughness of x. If we were interested in D2x(t), we might think of penalizing D4x(t) =4> cubic polynomials are not rough. • What about the weather data? We know it is periodic, and not very different from a sinusoid. • The Harmonic acceleration of x is Lx = uj2Dx+ D3x. with Lsin(cjx) = 0 = Lcos(cjx). • We can measure departures from a sinusoid by JL(x) = j [Lx(t)fdt. • Generally we can define m Lx{t) = Y,ak{t)Dkx{t). k=l Jan Koláček (SCI MUNI) M7777 Applied FDA Fall 2019 20 / 27 Smoothing Penalties Harmonic Acceleration for log A Jan Koláček (SCI MUNI) M7777 Applied FDA Fall 2019 21 / 27 Smoothing Penalties Harmonic Acceleration for log A Jan Koláček (SCI MUNI) M7777 Applied FDA Fall 2019 22 / 27 Smoothing Penalties Harmonic Acceleration for log A Jan Koláček (SCI MUNI) M7777 Applied FDA Fall 2019 23 / 27 Smoothing Penalties Choosing the Smoothing Parameter • Ordinary cross-validation ccvin)- 1 y^(r'-*(*'))2 W"A/^ (l-sn)2 ' 1=1 where Sa is the smoothing matrix. • Generalized cross-validation • GCV smooths more than OCV; even then, it may need to be tweaked a little to produce pleasing results. Jan Koláček (SCI MUNI) M7777 Applied FDA Fall 2019 24 / 27 Smoothing Penalties Both types of Cross-Validation for St. Johnes Precipitation Data CV for St. Johns -1 01 23456789 log lambda Jan Koláček (SCI MUNI) Fall 2019 25 / 27 Problems to solve O Melanoma Data • Load the variable melanoma from the f da package and plot it. • Fit the data using a B-spline basis of order 6 and a harmonic acceleration penalty. Try some values of A to optimize GCV. You will need to guess at the period to use; how does doubling and halving the period change the degrees of freedom at the optimal value of A? • Plot the velocity versus acceleration curves (phaseplanePlot command, Figure 1) for the fit using a Fourier basis and using the B-spline basis with a harmonic acceleration penalty. Are they substantially different? Do they provide evidence of sub-cycles? 0 Canadian Weather Data • Load the variable CanadianWeather from the f da package and select temperature data observed in Edmonton, Halifax, Montreal and Ottawa. • Fit temperature with Fourier bases and harmonic acceleration penalties at a number of values of A. • Plot GCV in terms of A and add the mean GCV. Choose the smooth for temperature that gives you the minimum mean gcv. Jan Koláček (SCI MUNI) M7777 Applied FDA Fall 2019 26 / 27 Problems to solve Jan Koláček (SCI MUNI) M7777 Applied FDA Fall 2019 27 / 27