A Semiparametric Nonlinear Mixed Effects Model with Penalized Splines Using Automatic Differentiation

Deze paper introduceert een schattingsprocedure voor semiparametrische niet-lineaire gemengde modellen met gepenaliseerde splines en automatische differentiatie, die in simulaties en een casestudie over de lengtegroei van zuigelingen een verbeterde inferentiële prestatie en een verlaagde computebelasting oplevert ten opzichte van bestaande methoden.

De kern

Stel je voor dat je een groep kinderen volgt die opgroeien. Je wilt weten hoe ze gemiddeld groeien (het "populatiepatroon"), maar je ziet ook dat elk kind zijn eigen tempo en timing heeft. Sommige kinderen groeien sneller, sommigen later, en sommigen zijn gewoon wat langer of korter.

Deze wetenschappelijke paper beschrijft een nieuwe, slimme manier om al deze verschillende groeicurves tegelijkertijd te analyseren. Hier is de uitleg in gewone taal, met een paar leuke vergelijkingen:

1. Het Probleem: De "Groeikunst"

Stel je voor dat je een schilderij wilt maken van hoe kinderen groeien.

• De oude manier: Je probeerde een vaste formule te vinden (bijvoorbeeld een rechte lijn of een simpele kromme) die bij iedereen past. Maar kinderen zijn complex; hun groei is geen rechte lijn.
• Het nieuwe idee: In plaats van een starre formule, gebruiken de auteurs een "slang" (een wiskundig iets dat ze penalized splines noemen). Deze slang is flexibel genoeg om de echte vorm van de groei na te bootsen, zonder dat je van tevoren weet hoe die eruit ziet.

2. De Uitdaging: Iedereen is anders

Elk kind volgt die "slang", maar op zijn eigen manier:

• Soms is de slang iets hoger of lager (lengteverschil).
• Soms is de slang versneld of vertraagd (snelheid van groei).
• Soms is de slang iets opgeschoven in de tijd (een kind dat later begint te groeien).

De auteurs noemen dit transformatieparameters. Het is alsof je de slang vasthoudt en hem uitrekt, verschuift of draait om hem op elk specifiek kind te laten passen.

3. De Oplossing: De "Slimme Slag" (De Methode)

De auteurs hebben een nieuwe rekenmethode bedacht om deze slang en de verschuivingen tegelijkertijd te berekenen. Ze gebruiken twee belangrijke trucjes:

A. De "Gladde Slang" (Penalized Splines)
Stel je voor dat je de slang te veel laat wiebelen; dan wordt hij onleesbaar. Je wilt dat hij glad is, maar wel de juiste vorm heeft.

• In het verleden moesten onderzoekers zelf kiezen hoe glad de slang moest zijn (een soort "gok").
• De nieuwe truc: De computer berekent automatisch hoe glad de slang moet zijn, net zoals hij de variatie tussen kinderen berekent. Het is alsof de computer zelf beslist hoeveel "stevigheid" de slang nodig heeft om niet te veel te wiebelen, maar wel de vorm te volgen.

B. De "Super-Rekenmachine" (Automatic Differentiation)
Het berekenen van deze modellen is als het proberen te vinden van het hoogste punt in een mistig berglandschap.

• Oude methode: Je loopt rond en voelt met je handen hoe het terrein eruitziet. Dat is traag en je kunt makkelijker een foutje maken.
• Nieuwe methode: Ze gebruiken een technologie genaamd Automatic Differentiation (via een softwaretool genaamd TMB). Dit is alsof je een drone hebt die direct en exact weet waar het hoogste punt is en hoe steil de helling is, zonder dat je zelf hoeft te rekenen. Het maakt de berekening veel sneller en nauwkeuriger.

4. Wat leverde dit op? (De Resultaten)

De auteurs hebben hun methode getest op twee manieren:

1. Simulaties (Virtuele kinderen): Ze maakten duizenden virtuele datasets. Hun nieuwe methode ("snmmTMB") was:

   • Sneller: Het deed er minder tijd over dan de oude methoden.
   • Nauwkeuriger: De "veiligheidsmarges" (betrouwbaarheidsintervallen) die ze berekenden, waren realistischer. De oude methode gaf soms te optimistische resultaten (te smalle marges), terwijl hun methode de onzekerheid eerlijker weergaf.

2. Echt Gebruik (Baby's): Ze pasten het toe op echte data van baby's die in hun eerste twee jaar gemeten werden.

   • Ze konden zien dat jongens gemiddeld iets langer zijn dan meisjes bij de geboorte.
   • Ze zagen dat als een baby te vroeg of te laat geboren is, de groeicurve zich bijna 1-op-1 verschuift in de tijd.
   • Ze kregen mooie, vloeiende grafieken die precies laten zien hoe baby's groeien, inclusief een "veiligheidsband" die laat zien hoe zeker we zijn van die lijn.

Samenvatting in één zin

Deze paper introduceert een snellere, slimmere manier om complexe groeicurve's te analyseren, waarbij de computer automatisch bepaalt hoe "glad" de lijnen moeten zijn en precies berekent hoe zeker we zijn van de resultaten, zonder dat de onderzoeker handmatig ingewikkelde formules hoeft te bedenken.

Het is alsof je van een handgetekende schets bent gegaan naar een 3D-scanner die de vorm perfect en snel vastlegt, inclusief een labeltje dat zegt: "Dit is 95% zeker de juiste vorm."

Technische samenvatting

Hier is een gedetailleerde technische samenvatting van het artikel "A Semiparametric Nonlinear Mixed Effects Model with Penalized Splines Using Automatic Differentiation" in het Nederlands.

Probleemstelling

Langeitudinale data (herhaalde metingen over tijd) komen veel voor in onderzoeksvelden zoals geneeskunde en psychologie. Een veelvoorkomend kenmerk is dat individuele trajecten een vergelijkbare onderliggende vorm hebben, maar verschillen in schaal, timing of andere subject-specifieke eigenschappen.
De auteurs richten zich op Semiparametrische Niet-lineaire Gemengde Effectenmodellen (SNMM), zoals voorgesteld door Ke en Wang (2001). In deze modellen wordt de populatietrajectorie geschat via een onbekende functie $f$, terwijl individuele afwijkingen worden gemodelleerd via transformatieparameters (random effects).

De huidige methoden voor het schatten van deze modellen (zoals de R-pakketten assist en sit) hebben echter belangrijke beperkingen:

1. Scheiding van schattingsprocedures: Bestaande algoritmen scheiden vaak de schatting van de vormfunctie van de schatting van de vaste en random effecten. Dit garandeert niet dat de procedure convergeert naar een maximum van de gezamenlijke likelihood, wat leidt tot onnauwkeurige inferentie en onvoldoende rekening houdend met onzekerheid.
2. Rekenkundige last: Het gebruik van smoothing splines vereist vaak een basisdimensie gelijk aan het aantal datapunten, wat computertijd kostbaar maakt.
3. Scheiding van gladheidsparameter: De gladheidsparameter moet vaak apart worden geselecteerd in plaats van gezamenlijk te worden geschat met andere variantiecomponenten, wat de schaalbaarheid beperkt.

Methodologie

De auteurs stellen een nieuwe schattingsprocedure voor die de volgende technische componenten combineert:

1. P-splines in Gemengd Model Representatie:

   • De populatietrajectorie $f$ wordt geschat met behulp van penalized splines (P-splines).
   • Door gebruik te maken van de eigendecompositie van de penalty-matrix, worden de spline-coëfficiënten opgesplitst in een deel dat als vaste effecten wordt behandeld (niet-gepenaliseerd) en een deel dat als random effecten wordt behandeld (gepenaliseerd).
   • Dit stelt de auteurs in staat om de gladheidsparameter ($\lambda$) direct te schatten als een variantiecomponent binnen het gemengde model, in plaats van deze apart te selecteren.

2. Laplace-benadering voor Marginalisatie:

   • De likelihood vereist integratie over de random effects (zowel de subject-specifieke transformaties als de spline-random effects). Omdat deze integratie analytisch niet oplosbaar is vanwege de niet-lineariteit, wordt de Laplace-benadering gebruikt.
   • Dit benadert de integraal door een lokale Gaussische benadering rond het maximum van de integrand.

3. Automatische Differentiatie (AD):

   • Een cruciale innovatie is het gebruik van Template Model Builder (TMB) in R. TMB maakt gebruik van C++ en operator-overloading om automatische differentiatie toe te passen.
   • Hierdoor worden exacte eerste en tweede orde afgeleiden (gradienten en Hessiaan) van de likelihood-functie gegenereerd met machine-precisie, zonder dat deze handmatig moeten worden afgeleid. Dit maakt de methode flexibel voor complexe functies en versnelt de convergentie.

4. Knooppunten (Knots) Selectie:

   • Omdat de transformatieparameters de x-as beïnvloeden, is de locatie van de knooppunten niet direct observeerbaar. De auteurs stellen een procedure voor waarbij de knooppunten worden gefixeerd op een geschaald domein (bijv. [0,1]) dat dynamisch wordt aangepast aan de geschatte variantie van de random effects, zodat differentieerbaarheid behouden blijft.

Belangrijkste Bijdragen

• Gezamenlijke Schatting: De methode schat de gladheidsparameter, de vaste effecten en de variantiecomponenten van de random effects simultaan via één gezamenlijke likelihood.
• Efficiëntie en Nauwkeurigheid: Door AD en TMB te gebruiken, wordt de berekening van afgeleiden en de optimalisatie aanzienlijk versneld en nauwkeuriger dan bij traditionele methoden.
• Flexibiliteit: De aanpak is implementatie-vriendelijk en kan worden uitgebreid met andere verdelingen of covariantiestructuren zonder de onderliggende differentiatie-logica te hoeven herschrijven.
• Robuuste Inferentie: Het biedt een kader voor het construeren van betrouwbaarheidsintervallen voor zowel de populatiecurve als de individuele trajecten, inclusief simultane banden.

Resultaten

De methode is getest in twee simulatiestudies en toegepast op een reëel dataset:

1. Simulatiestudies:

   • De nieuwe methode (snmmTMB) werd vergeleken met het bestaande pakket assist.
   • Resultaat: snmmTMB leverde een betere dekking (coverage) van de betrouwbaarheidsbanden op, dichter bij de nominale waarde (bijv. 95%), zelfs bij hoge variantie en kleine steekproeven.
   • De betrouwbaarheidsbanden van snmmTMB waren smaller en stabieler dan die van assist.
   • Rekentijd: snmmTMB was aanzienlijk sneller (bijv. 5-40 seconden versus 7-170 seconden voor assist) en toonde minder variabiliteit in rekentijd.

2. Toepassing op Infanthoogte (SMOCC-data):

   • De methode werd toegepast op lengtemetingen van Nederlandse kinderen tot twee jaar.
   • Het model schatte een gladde populatiegroei-curve met covariaten voor geslacht en zwangerschapsduur.
   • De resultaten waren consistent met bekende groeipatronen (snelle groei in de eerste 6 maanden).
   • Een parametrische bootstrap analyse bevestigde dat de Laplace-benadering adequaat was voor deze data, aangezien de asymptotische standaardfouten goed overeenkwamen met de empirische variantie uit de bootstrap.

Betekenis en Conclusie

Dit artikel biedt een oplossing voor een langdurig probleem in de statistische modellering van niet-lineaire gemengde effectenmodellen. Door de combinatie van P-splines, de gemengde model-representatie en automatische differentiatie via TMB, overwinnen de auteurs de beperkingen van eerdere methoden.

De belangrijkste implicaties zijn:

• Betrouwbaardere inferentie: Onderzoekers kunnen nu meer vertrouwen hebben in de schattingen van populatietrajecten en hun onzekerheid.
• Scalabiliteit: De methende is geschikt voor grotere datasets dankzij de efficiëntie van AD.
• Praktische toepasbaarheid: De open-source implementatie (beschikbaar via GitHub) maakt de methode direct toepasbaar voor onderzoekers in biostatistiek, psychologie en andere domeinen met complexe longitudinale data.

De auteurs sluiten af met de opmerking dat toekomstig werk gericht kan zijn op het verbeteren van modelselectie voor random effects (bijv. via aangepaste AIC) en het toepassen van vormbeperkingen (zoals monotonie) op de spline-curves.