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 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 10.0 7.5 o -i—' CO -i—i 'Cl ö 5.0 CD l_ Q_ 2.5 tM Ail 100 200 Days 300 Jan Koláček (SCI MUNI) M7777 Applied FDA Fall 2019 10 / 27 Smoothing Penalties Smoothing with In A St. Johns • ( » • • • m • • * -• • • • F vi • • • •• • • m ■iT. • • • • •A* t7tJ Mf*W\/ !•./< » 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 • » • m • • mm • A • • • • • If Vhf», • • • • m ,* • • • - • • • • 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 • ( » • M • • • • • • • • • * • • • If Vhr», • • • • • ■ m - • • - • • • • 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 • •• • • • m • • • A • ( I • • • • • •• • • • • • A« • • - • • • • • ► *S « • • 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 • ( » • • • M • • • » • • • m ( 1 • • » • m • • • •• • * ■ • j • -V « • ^AA« Vi • • • • • • • • * * * • tl • ( ) 1( • )0 2( )0 3( )0 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) + c;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 4 3 7 l'l 15 ~m 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