A review of spline function procedures in R

Termin 'spline’ odnosi się do narzędzia rzemieślniczego, giętkiego cienkiego paska drewna lub metalu, używanego do kreślenia gładkich krzywych. Kilka ciężarków byłoby przyłożonych w różnych pozycjach, więc pasek zginałby się zgodnie z ich liczbą i pozycją. W ten sposób pasek był zmuszany do przejścia przez zestaw stałych punktów: metalowe kołki, żebra łodzi itp. Na płaskiej powierzchni były to często ciężarki z przymocowanym hakiem i w ten sposób łatwe do manipulowania. Kształt wygiętego materiału w naturalny sposób przyjąłby formę krzywej splajnu. Podobnie, splajny są używane w statystyce, aby matematycznie odtworzyć elastyczne kształty. Węzły są umieszczone w kilku miejscach w zakresie danych, aby zidentyfikować punkty, w których sąsiednie elementy funkcjonalne łączą się ze sobą. Zamiast metalowych lub drewnianych pasków, do dopasowania danych pomiędzy dwoma kolejnymi węzłami wybierane są gładkie elementy funkcyjne (zwykle wielomiany niskiego rzędu). Typ wielomianu oraz liczba i umiejscowienie węzłów jest tym, co następnie definiuje typ splajnu.

Przykład motywujący

Po wprowadzeniu uogólnionych modeli addytywnych (GAM) w 1986 roku, zastosowanie modelowania splajnu stało się uznanym narzędziem w statystycznej analizie regresji. Aby to zilustrować, rozważ dane na temat zestawu 892 kobiet poniżej 50 lat zebranych w trzech wioskach w Afryce Zachodniej (dane dostępne w pliku dodatkowym 1: dodatek). Chcielibyśmy zbadać związek między wiekiem (w latach) i surowej miary tkanki tłuszczowej, która jest grubość fałdu skórnego triceps. Rycina 1 pokazuje zależność między wiekiem a grubością triceps skinfold mierzoną w skali logarytmicznej. Więcej informacji na temat danych można znaleźć na stronie .

Prosty model regresji w postaci yi=β0+β1xi+ε,i=1,…,n, z trudem dałby przybliżenie obserwowanego wzoru, ponieważ oczywiste jest, że zależność nie jest liniowa. Model można rozszerzyć, aby uwzględnić efekty nieliniowe za pomocą pewnych wielomianów. Wówczas efekty nieliniowe mogłyby być modelowane wielomianem stopnia 3 danym przez:

$$ y_{i}=}alfa_{0}+alfa_{1} u_{i}+alfa_{2} u_{i}^{2}+alfa_{3} u_{i}^{3}+epsilon $$

(1)

gdzie u jest funkcją x zwaną funkcją bazową, zdefiniowaną tutaj przez:

$$U= funkcja bazowa $$

Model regresji opisany równaniem 1 jest nadal modelem liniowym, pomimo tego, że podaje nieliniową funkcję zmiennej predykcyjnej. Model ten jest nadal liniowy w zakresie współczynników i może być dopasowany przy użyciu zwykłych metod najmniejszych kwadratów. Podstawę można utworzyć w R za pomocą funkcji poly(x,3) z danymi wejściowymi x (odnoszącymi się do zmiennej), oraz p (odnoszącymi się do stopnia wielomianu). Prowadzi to do prostego jednowariantowego gładkiego modelu postaci: yi=f(xi)+ε gdzie f() jest pewną funkcją/transformacją predyktora. Taki model można łatwo dopasować w programie R za pomocą: lm(y ∼poly(x,3)). Pomimo prostoty, regresja wielomianowa ma kilka wad, z których najważniejszą jest nielokalność. Oznacza to, że dopasowana funkcja przy danej wartości x0 zależy od wartości danych odległych od tego punktu. Można to łatwo zobaczyć w działaniu, dopasowując wielomian do zbioru danych i przesuwając jeden z punktów danych w pobliżu prawej krawędzi w górę lub w dół. W rezultacie dopasowana funkcja zwykle zmieni się daleko od tej współrzędnej x.

Rozważmy, zamiast dopasowania globalnego wielomianu, podział zakresu x na mniejsze przedziały, wykorzystując arbitralną liczbę i położenie punktów, τ, zwanych również węzłami. Prosty model ciągły cząstkowy można dopasować definiując funkcje: f1(x)=1,f2(x)=x,f3(x)=(x-τ1)+,f4(x)=(x-τ2)+,…, przy czym „+” jest funkcją zdefiniowaną jako: $$

Zbiór tych funkcji prowadzi do funkcji zespolonej f(x).

Definicja splajnów

Metalowy splajn rysownika może przyjmować dowolne kształty, na przykład przekrój skrzydła samolotu lub spirali pompy odśrodkowej. W zastosowaniach statystycznych przyjmiemy krzywe w postaci f(X), tj. jedną wartość y dla każdego x. Predyktor x może być pojedynczą zmienną lub wieloma zmiennymi. W naszej dyskusji skupimy się prawie wyłącznie na funkcji jednozmiennowej o \(X w {R}). Zdefiniuj zbiór węzłów τ1<…<τK w przedziale X. Splajnem f(X) będzie funkcja gładka, spełniająca pewne własności różniczkowalności, o których będzie mowa poniżej, taka, że f(X) jest wielomianem stopnia d. Splajny drewniane lub metalowe mają ciągłe pochodne wszystkich rzędów, ponieważ są obiektem fizycznym. Nie jest to prawdziwe dla splajnów statystycznych. Narzucamy raczej kryterium gładkości, że wszystkie pochodne rzędu mniejszego niż d są ciągłe. Splajn fizyczny jest liniowy poza ostatnim węzłem i możemy nałożyć dodatkowe ograniczenie, że pochodne rzędu 2 lub większego są równe zero w najbardziej lewym i najbardziej prawym węźle; splajny z tym dodatkowym ograniczeniem są znane jako splajny „ograniczone” lub „naturalne”. W celu uzyskania bardziej elastycznych krzywych można zwiększyć liczbę węzłów lub stopień wielomianu. Istnieje jednak kompromis; zwiększenie liczby węzłów może spowodować niedopasowanie do danych i zwiększyć wariancję, podczas gdy zmniejszenie liczby węzłów może skutkować sztywną i restrykcyjną funkcją, która ma większą tendencyjność.

Przedstawienie za pomocą funkcji bazowych

Załóżmy, że nieznana funkcja f jest reprezentowana przez funkcję splajnu o stałej sekwencji węzłów i stałym stopniu d. Ponieważ te ostatnie funkcje tworzą przestrzeń wektorową V, można zapisać f jako

$$ f(X)=

B_{k} (X), $$

(2)

gdzie Bk to zbiór funkcji bazowych definiujących V, a βk to związane z nimi współczynniki spline. Przy węzłach istnieje k+1 wielomianów stopnia d wraz z d∗k ograniczeniami, co prowadzi do (d+1)(k+1)-d∗k=d+k+1 wolnych parametrów; dla naturalnego splajnu istnieje k wolnych parametrów. Ponieważ βB=(βA)(A-1B)=γB∗ dla dowolnej niesingularnej macierzy A istnieje nieskończenie wiele możliwych zbiorów bazowych dla dopasowania splajnu.

Przedstawienie w (2) ma tę zaletę, że estymacja f sprowadza się do estymacji współczynników βk. Dokładniej, wyrażenie w (2) jest liniowe w wektorze współczynników β=(β1,…,βK+d+1). Zatem estymacja f może być postrzegana jako problem optymalizacyjny, który jest liniowy w przekształconych zmiennych B1(X),…,BK+d+1(X), co pozwala na wykorzystanie dobrze ugruntowanych technik estymacji do użycia splajnów w szerokim zakresie (uogólnionych) modeli regresji wielozmiennowej. Co ważne, modelowanie splajnowe redukuje estymację funkcji f() do estymacji niewielkiego zbioru współczynników o rzeczywistej wartości.

Jak wskazują różni autorzy (np. wysoka elastyczność modelowania splajnowego ma swoją cenę w postaci szeregu parametrów dostrajania. Dwa z nich, wybór funkcji bazowych B i stopień d wielomianów bazowych, okazują się mieć niewielki wpływ. Wielomiany sześcienne (d=3) są zwykle standardem, ponieważ dają one krzywe, które wydają się idealnie gładkie dla ludzkiego oka. Jeśli pochodne dopasowanych krzywych są interesujące, wyższy rząd jest czasami odpowiedni, ale ogólnie rzecz biorąc, dopasowania dla d>3 są efektywnie nierozróżnialne. Dopasowania z d=1 lub d=2 mają prawie identyczne właściwości statystyczne, ale będą wydawać się bardziej poszarpane. Wybór pomiędzy dwoma zestawami podstaw B i B∗ z definicji nie zmieni przewidywań z dopasowania, a więc sprowadza się do kwestii wygody.

Dwa kluczowe wybory dotyczą liczby i rozstawu węzłów oraz użycia (lub nie) funkcji kary, np. zintegrowanej drugiej pochodnej splajnu. Gdy nie ma kary, tworzenie zmiennych przekształconych może być wykonane oddzielnie, a nowe zmienne są po prostu włączone do standardowego dopasowania modelu; nie jest wymagana żadna modyfikacja podstawowej procedury regresji. Podejście to jest często określane jako splajny regresyjne; elastyczność wynikowej funkcji nieliniowej jest całkowicie funkcją liczby węzłów. Z drugiej strony, włączenie kary za wygładzanie wymaga modyfikacji procedury dopasowania w celu jej uwzględnienia. Musi ona być uwzględniona w każdej funkcji regresji osobno. Wynikowe splajny wygładzające mają kilka pożądanych właściwości, ale dodatkowa złożoność funkcji wygładzającej może być powodem, dla którego nie są częściej wykorzystywane w zastosowaniach.

Chociaż przeprowadzono znaczne badania w celu zbadania właściwości matematycznych różnych podejść do splajnów (patrz , statystycy stosowani i analitycy danych wydają się prawie nie być świadomi tych wyników podczas korzystania z modelowania splajnów w praktycznych zastosowaniach. W rzeczywistości, wiele artykułów zidentyfikowanych przez naszą kwerendę internetową nie zawierało uzasadnienia wyboru zastosowanej metody spline.

Popularne podstawy spline

Istnieje wiele opcji definicji funkcji bazowych Bk, gdzie różne podstawy spline różnią się pod względem ich właściwości numerycznych . W tym rozdziale przedstawimy niektóre z najbardziej popularnych podstaw splajnu, a mianowicie podstawę obciętego szeregu potęgowego, podstawę splajnu B oraz podstawę splajnu kardynalnego.

Ocięty szereg potęgowy i splajny sześcienne

Podstawa obciętego szeregu potęgowego jest zdefiniowana przez funkcje bazowe

$$B_{1}(x) = 1, B_{2}(x) = x,……, B_{d+1}(x) = x^{d}, $$

$$B_{d+2}(x) = (x- -tau_{1})_{+}^{d},…, B_{K+d+1} = (x -tau_{k})_{+}^{d} $$

Zaletą powyższych funkcji bazowych jest ich łatwa interpretacja: Zaczynając od „podstawowego” wielomianu stopnia d określonego na (pierwszy wiersz równania), odchylenia od wielomianu podstawowego są kolejno dodawane do funkcji splajnu na prawo od każdego z K węzłów (drugi wiersz). Splajn o obciętej podstawie potęgowej jest d-1 razy różniczkowalny w węzłach i ma d+K stopni swobody. Użytkownik może stosunkowo łatwo utworzyć obcięty szereg potęgowy w R. Niech x reprezentuje pewne obserwacje w , wtedy obcięta podstawa potęgowa stopnia d=3 z 5 węzłami równo rozmieszczonymi wzdłuż zakresu x może być utworzona przy użyciu kodu 1 w pliku dodatkowym 1: Dodatek (Rys. 2).

Cechą obciętych szeregów potęgowych jest to, że podpory funkcji nie są lokalne, a niektóre z Bk są określone w całym zakresie danych . Może to prowadzić do wysokich korelacji pomiędzy niektórymi splajnami bazowymi, implikując niestabilności numeryczne w estymacji splajnu. Dla obciętej podstawy szeregu potęgowego, przykład jest podany w , Rozdział 5.

Splajny sześcienne są tworzone przez użycie wielomianu sześciennego w przedziale pomiędzy dwoma kolejnymi węzłami. Splajn ma cztery parametry na każdym z K+1 regionów minus trzy ograniczenia dla każdego węzła, co daje K+4 stopnie swobody.

Funkcja splajnu sześciennego, z trzema węzłami (τ1,τ2,τ3) będzie miała 7 stopni swobody. Korzystając z reprezentacji podanej w Eq. 2, funkcję można zapisać jako:

$$ f(X)= \beta_{0} + \beta_{1} X + \beta_{2} X^{2} + \beta_{3} X^{3} + \beta_{4} (X-tau_{1})^{3} + \beta_{5} (X-tau_{2})^{3} + \beta_{6} (X-tau_{2})^{3} + \beta_{6} (X-tau_{3})^{3} $$

B-splajny

Podstawa B-splajnu jest powszechnie stosowaną podstawą splajnu, która jest oparta na specjalnej parametryzacji splajnu sześciennego. Podstawa B-splajnu , jest oparta na sekwencji węzłów

$$begin{aligned} \xi_{1} ∗ ∗ ∗ ∗ ∗ ∗ &le \u_xi_{d+1} < \xi_{d+2} < \u002600↩ \u002600↩ \u002600↩ \u002600↩ \u002600↩ \u002600↩ \u002600↩ \\ ^8799>< ^i_{d + K + 2} ^i_{d + K + 3} ^i_{d + K + 2} ^i_{d + K + 2} ^i_{d + K + 2} ^,, ^end{aligned} $$

gdzie zbiory ξd+2 := τ1,…,ξd+K+1:=τK oraz ξd+1:=a,ξd+K+2:=b nazywane są odpowiednio „węzłami wewnętrznymi” i „węzłami brzegowymi”. Wybór dodatkowych węzłów ξ1,…,ξd i ξd+K+3,…,ξ2d+K+2 jest w zasadzie dowolny. Powszechnie przyjętą strategią jest ustawienie ich na równi z węzłami brzegowymi. Alternatywnie, jeśli węzły wewnętrzne i brzegowe ξd+1<…<ξd+K+2 są tak dobrane, by były równomiernie oddalone, tzn, ξk+1-ξk=δ ∀k∈{d+1,…,d+K+1}, węzły graniczne można umieścić w punktach ξd+1-δ,…,ξd+1-d-δ i ξd+K+2+δ,…,ξd+K+2+d-δ.

Dla d>0, funkcje bazowe B-spline stopnia d (oznaczane przez $B_{k}^{d}(x)$ są określone wzorem rekurencyjnymPrzypis 1

$$ \begin{aligned} B_{k}^{d}(x)&=\frac{x-\xi_{k}}{\xi_{k+d}-\xi_{k}}B_{k}^{d-1}(x)-\frac{\xi_{k+d+1}-x}{\xi_{k+d+1}-\xi_{k+1}}B_{k+1}^{d-1}(x),\\k &= 1,….,K+d+1, ^end{aligned} $$

gdzie

$$B_{k}^{0}(x)= ^begin{array}{cc} 1, & ^xi_{k} \\ x < \xi_{k+1} \ 0, & \tekst{else} \ end{array} \prawda. $$

oraz $B_{k}^{0}(x) $, jeśli ξk=ξk+1. B-splajny mają tę zaletę, że funkcje bazowe mają lokalne wsparcie. Mówiąc dokładniej, są one większe od zera w przedziałach rozpiętych przez d+2 węzły i zerowe gdzie indziej. Ta własność skutkuje wysoką stabilnością numeryczną, a także efektywnym algorytmem konstrukcji funkcji bazowych, zobacz szczegóły.

Naturalne splajny sześcienne i kardynalne

Splajny wielomianowe, takie jak sześcienny lub B-splajn, mogą być nieregularne na granicach danych. Aby rozwiązać ten problem, splajny naturalne są splajnami sześciennymi, które mają dodatkowe ograniczenie, że są liniowe w ogonach węzłów granicznych (-∞,a], Rozdział 4.

Oprócz obciętych szeregów potęgowych splajnów naturalnych, B-splajnów i podstaw splajnów kardynalnych istnieją różne inne – mniej popularne – podstawy. W celu zapoznania się z nimi odsyłamy do książek autorstwa .

Penalizowane splajny

Przedstawione do tej pory splajny są często określane jako splajny regresyjne. Oprócz wyboru podstawy splajnu (B-splajn, obcięty szereg potęgowy, itp.), należy wybrać liczbę węzłów i położenie węzłów. Oczywiście, te parametry dostrajania mogą mieć istotny wpływ na estymowany kształt funkcji splajnu: Duża liczba węzłów oznacza dużą elastyczność, ale może również prowadzić do nadmiernego dopasowania do danych. I odwrotnie, mała liczba węzłów może skutkować „nadmiernie wygładzoną” estymatą, która jest podatna na niedopasowanie (patrz ).

Popularnym podejściem ułatwiającym wybór pozycji węzłów w modelowaniu splajnu jest użycie splajnów karanych. Biorąc pod uwagę i.i.d. próbkę danych (x1,y1),…(xn,yn), splajn penalizowany jest rozwiązaniem problemu

$$hat{beta} = \tekst{argmax}_{beta} \left \,, $$

gdzie lβ oznacza log-likelihood (lub, w przypadku regresji Coxa, częściowe log-likelihood), a Jr jest karą za szorstkość, która staje się mała, jeśli funkcja splajnu jest „gładka”. Generalnie, penalizowane splajny opierają się na założeniu, że nieznana funkcja f jest modelowana przez splajn o dużej liczbie węzłów, co pozwala na wysoki stopień elastyczności. Z drugiej strony, przybliżona estymata splajnu, która ma wysoką wartość lβ i jest bliska wartościom danych, skutkuje dużą wartością Jβ. Maksymalizacja tej funkcji implikuje zatem kompromis pomiędzy gładkością a dopasowaniem modelu, który jest kontrolowany przez parametr dostrajania λ≥0.

Specjalnym przypadkiem jest penalizowany problem najmniejszych kwadratów

$$ \hat{\beta} = \text{argmin}_{\beta} \left $$

(3)

w regresji gaussowskiej. Kara $J_{beta } \, \int _{a}^{b} \left (\partial ^{2} f / \partial x^{2}right)^{2} dx $ wyraża „gładkość” funkcji splajnu pod względem drugiej pochodnej f. Dla danego λ można pokazać, że rozwiązaniem jest naturalny splajn sześcienny z sekwencją węzłów x(1)<…<x(n), tzn, pozycje węzłów nie muszą być wybierane, lecz są „naturalnie” dane przez uporządkowane unikalne wartości danych X. W literaturze ten typ splajnu jest określany jako smoothing spline . Można wykazać, że splajn wygładzający interpoluje dane, jeśli λ=0, podczas gdy λ=∞ implikuje funkcję liniową. Zauważmy, że splajny wygładzające są szczególnym przypadkiem bardziej ogólnej klasy cienkich splajnów płytowych , które pozwalają na rozszerzenie kryterium w równaniu (3) na wyższe wymiary xi (patrz , Rozdział 4.15], i dla szczegółów).

Wygodną własnością splajnów wygładzających jest to, że kara Jβ może być zapisana jako β⊤Ωβ z odpowiednio zdefiniowaną macierzą kary Ω. Zatem rozwiązanie (3) jest dane przez estymatę karną metodą najmniejszych kwadratów

$$ ^hat{beta} = ^left(B^{top} B + ^lambda ^Omega}right)^{-1} B^{{top} y $$

(4)

gdzie B jest macierzą o wymiarze n×n zawierającą funkcje bazowe splajnu naturalnego obliczone dla wartości danych. Wektor y zawiera wartości odpowiedzi y1,…,yn. W praktyce istnieją bardzo wydajne algorytmy do obliczania ∗(∗beta) w (4) . Zamiast określać naturalną podstawę splajnu dla f, możliwa jest praca z nieograniczoną podstawą B-splajnu, ponieważ kara w (3) automatycznie nakłada ograniczenia liniowości w węzłach x(1) i x(n) (patrz , Rozdział 5, i , Rozdział 2). Jeśli chodzi o podstawę B-spline, wyniki estymacji nie będą zależały od wyboru węzłów brzegowych: możliwe jest albo użycie x(1) i x(n) jako węzłów brzegowych, albo włączenie x(1) i x(n) do zbioru węzłów wewnętrznych.

Jeśli n jest duże, a przedział jest gęsto pokryty obserwowanymi danymi, zwykle nie jest konieczne umieszczanie węzła przy każdym xi,i=1,…,n. Zamiast tego, splajn wygładzający może być aproksymowany przez penalizowany splajn regresyjny, który wykorzystuje zredukowany zbiór węzłów. Bardzo popularną klasą penalizowanych splajnów regresyjnych są P-splajny, które bazują na sześciennej podstawie B-splajnu i na „dużym” zbiorze równomiernie oddalonych węzłów (zwykle 10-40). Zamiast obliczać całkę w (3), P-spliny są oparte na karze różnicy drugiego rzędu zdefiniowanej przez

$$J^{*}_{beta} = ^suma_limits_{k=3}^{K+4} \left(\Delta^{2} \beta_{k} \right)^{2} \, $$

co w przypadku równomiernie rozmieszczonych węzłów można wykazać jako przybliżenie do Jβ. Operator różnicy drugiego rzędu Δ2 jest zdefiniowany przez Δ2βk:=(βk-βk-1)-(βk-1-βk-2). Karę można zatem wyrazić jako β⊤Pβ, gdzie P jest zdefiniowane przez D⊤D, a D jest macierzą różnic. Można łatwo wywnioskować, że wynikowy estymator β ma taką samą strukturę jak 2, z Ω zastąpionym przez P.

Wygodną własnością P-splajnów jest to, że są one stabilne numerycznie oraz bardzo łatwe do zdefiniowania i implementacji. W szczególności, znacznie łatwiej jest utworzyć macierz różnicy D niż macierz Ω. Ponadto, łatwo jest rozszerzyć karę Jβ (a tym samym macierz D) na różnice wyższego rzędu Δq z q>2. Możliwe jest również użycie sekwencji węzłów, które nie są równomiernie rozłożone; w tym przypadku należy wprowadzić wagi. Ponieważ P-splajny z nierównomiernie rozmieszczonymi węzłami są rzadko stosowane w praktyce, nie rozważamy ich tutaj i odwołujemy się do innych rozwiązań.

Smoothing splines i P-splajny do pewnego stopnia przezwyciężają problem wyboru węzłów. Ich filozofia polega na użyciu dużej liczby węzłów, a następnie pozwoleniu, aby λ kontrolowała ilość gładkości. Powoduje to powstanie jednego dodatkowego parametru do dostrajania, przy czym nie ma ogólnej zgody co do sposobu dostrajania tego parametru. Niektóre popularne sposoby określania „optymalnej” wartości λ wykorzystują uogólnioną walidację krzyżową (GCV), AIC lub reprezentację mieszanego modelu .

Splajny w R

Podstawowy pakiet instalacyjny R zawiera zestaw funkcji, które mogą dopasować proste wielomianowe splajny i splajny wygładzające. Dalsze funkcje są zawarte w bibliotece splines napisanej przez DM Batesa i WN Venablesa. Pakiet ten był przez wiele lat koniem roboczym dopasowywania splajnów i jest obecnie częścią podstawowej dystrybucji R. Istnieje ponad 100 innych pakietów, które podczas ładowania zależą od splajnów. Pakiet zawiera kilka funkcji do tworzenia podstaw splajnów, takich jak bs dla B-splajnów i ns dla naturalnych splajnów, które są powszechnie używane, ale także kilka bardziej wyspecjalizowanych funkcji do tworzenia funkcji bazowych (takich jak periodicSpline, która tworzy okresowe interpolacje splajnów) lub poleceń, które są przydatne, takie jak polecenie predict.Domyślne wartości bs utworzą sześcienną podstawę B-splajnu z dwoma węzłami brzegowymi i jednym węzłem wewnętrznym umieszczonym na medianie obserwowanych wartości danych. Użytkownik może uzyskać większą elastyczność poprzez zwiększenie rozmieszczenia i liczby węzłów i/lub zmianę ich lokalizacji. Rysunek 3 (kod 2 w pliku dodatkowym 1: Appendix) przedstawia B-splajny utworzone przy użyciu różnych opcji. Górna część przedstawia splajny liniowe, tj. wielomiany pierwszego rzędu (stopień jeden) połączone ze sobą na jednakowo odległych węzłach. Dolna część przedstawia wielomiany sześcienne (stopień 3).

Należy zauważyć, że B-splajny utworzone w R za pomocą bs() są automatycznie ograniczone zakresem danych, a dodatkowe węzły (τ1,…,τd) są ustawione na równi z węzłami brzegowymi, co daje wiele węzłów na obu końcach dziedziny. Takie podejście jest użyteczne w przypadkach jednowariantowych i ma pewne atrakcyjne obliczeniowo cechy. Jednakże, jeśli ktoś pracuje nad dwuwymiarowym problemem wygładzania, używając iloczynów tensorowych B-splajnów, lub gdy pracuje z P-splajnami, ta podstawa jest nieodpowiednia i może prowadzić do błędnych wyników.

Naturalne splajny mogą być tworzone w pakiecie splines, przy użyciu polecenia ns. Domyślnie, o ile użytkownik nie określi stopni swobody lub węzłów, funkcja zwraca linię prostą w obrębie węzłów brzegowych. Rysunek 4 (kod 3 w pliku dodatkowym 1: Appendix pokazuje naturalne splajny utworzone przy użyciu różnych opcji.

Aby zilustrować, jak te funkcje mogą być wykorzystane w praktyce, rozważmy ponownie dane z rozdziału 2.0.1. Rysunek 5 (utworzony przez (kod 4 w pliku dodatkowym 1: Appendix)) pokazuje dopasowania uzyskane przy użyciu następujących poleceń: poly() dla prostych ortogonalnych splajnów wielomianowych, smooth.spline() dla splajnów wygładzających, bs() i ns() z biblioteki splines, odpowiednio dla B-splajnów i splajnów naturalnych. Górny lewy wykres przedstawia proste liniowe dopasowanie do danych (linia przerywana) oraz dopasowanie wielomianowe trzeciego stopnia, które jest w stanie uchwycić bardziej złożone zależności pomiędzy zmiennymi. Wykres w prawym górnym rogu jest szczególnie interesujący, ponieważ przedstawia dopasowania przy użyciu domyślnych wartości funkcji splajnu. Zielona linia pochodzi z funkcji poly() oraz ns(), które domyślnie definiują linię prostą. Na drugim biegunie, niebieska linia jest dopasowaniem z funkcji smooth.spline(), która jeśli nie określono stopni swobody ma tendencję do niedogładzania danych, tj. produkuje bardzo elastyczne dopasowanie wiggly oparte – tutaj- na 45 stopniach swobody. Wizualnie rozsądne dopasowanie do danych można uzyskać po określeniu czterech stopni swobody (lewy dolny wykres). Można zauważyć, że istnieją pewne różnice w zależności od wybranej bazy. Podstawa wielomianowa (czarna linia) jest nieco bardziej elastyczna niż pozostałe, szczególnie w wyższych grupach wiekowych. Z drugiej strony, splajn wygładzający ograniczony do zaledwie czterech stopni swobody jest bardziej sztywny niż inne podejścia, ale prawdopodobnie nadmiernie wygładza dane w niskim wieku, między 0 a 10 rokiem życia. Pomiędzy tymi dwoma ekstremami, B-splajny i naturalne splajny zapewniają bardzo podobne dopasowania, które wychwytują efekt małego wieku i mają tendencję do pozostawania pod mniejszym wpływem skrajnych przypadków na końcu spektrum wiekowego. Na koniec, prawy dolny wykres pokazuje, jak bardzo elastyczne stają się dopasowania z dodatkowymi stopniami swobody i sugeruje potencjalne niedopasowanie z powodu użycia zbyt wielu stopni swobody.

A note on degrees of freedom

W praktyce zawsze warto zdefiniować splajn za pomocą stopni swobody. Takie podejście jest szczególnie przydatne podczas pracy z B-splajnami i splajnami naturalnymi. B-splajny mają d+K, podczas gdy funkcja bazowa naturalnego splajnu sześciennego z K węzłami ma odpowiednio K+1 stopni swobody. Domyślnie funkcja bs w R tworzy B-splajny stopnia 3 bez węzłów wewnętrznych i z węzłami granicznymi zdefiniowanymi na zakresie zmiennej X. Jako taka funkcja tworzy trzy funkcje bazowe. Rozważmy teraz następujący przypadek: gdy użytkownik zdefiniuje B-splajnę z węzłem wewnętrznym na medianie X (bs(x,knots=median(x))), program utworzy cztery funkcje (d=3 plus K=1 węzłów wewnętrznych, cztery stopnie swobody). Jeśli jednak użytkownik określi w funkcji węzły brzegowe w ramach argumentu knots (bs(x,knots=c(min(x),median(x),max(x))))), to funkcja będzie miała sześć stopni swobody (d=3 plus k=3). Podobną ostrożność należy zachować w przypadku funkcji ns.

Przy pracy ze splajnami wygładzającymi nie jest łatwo określić stopnie swobody, ponieważ będą się one zmieniać w zależności od wielkości kary. Jednak w praktyce penalizowane splajny mogą być również ograniczone do maksymalnej liczby stopni swobody lub pożądanych stopni swobody.

Inne pakiety splajnów

Szerzej mówiąc, rozszerzona lista pakietów splajnów zawiera albo podejścia dość podobne do przedstawionych tutaj, albo bardzo wyspecjalizowane przypadki, które są ukierunkowane na konkretne zastosowania. W Tabeli 1 przedstawiono niektóre z tych pakietów wraz z liczbą pobrań. Liczba ta odnosi się do liczby pobrań pakietu, ale nie do unikalnych użytkowników. Szczegółowy opis wszystkich tych podejść wykracza poza zakres tej pracy.

Tabela 1 Pakiety R używane do tworzenia splajnów

Arquidia Mantina