I detta dokument förespråkas en utvidgning av den roll som Bayesian statistiska slutsatser spelar när det formellt kvantifierar osäkerheten i datormodeller som definieras av system med vanliga eller partiella differentialekvationer. Vi antar perspektivet att implicit definierade oändliga dimensionella funktioner som representerar modelltillstånd är objekt som kan härledas probabilistiskt. Vi utvecklar en allmän metod för probabilistisk integration av differentialekvationer genom modellbaserad uppdatering av en gemensam föregående åtgärd om funktionernas utrymme och deras tidsbundna och rumsliga derivat. Detta resulterar i ett bakre mått över funktioner som återspeglar hur väl de uppfyller systemet med differentialekvationer och motsvarande initial- och gränsvärden. Vi visar hur denna bakre åtgärd naturligt kan införlivas i Kennedys och O'Hagans ramverk för osäkerhetskvantifiering och ger ett helt Bayesianskt förhållningssätt till modellkalibrering. Genom att ta denna probabilistiska ståndpunkt, kan den fulla kraften av Bayesian inference utnyttjas när man försöker konsekvent kvantifiera och sprida epistemisk osäkerhet i datormodeller av komplexa naturliga och fysiska system. Ett brett utbud av exempel ges för att illustrera potentialen i detta ramverk för att karakterisera diskretiseringsosäkerhet, inklusive initialvärde, fördröjning och gränsvärde differentialekvationer, samt partiella differentialekvationer. Vi demonstrerar också vår metodik på ett storskaligt system, genom att modellera diskretization osäkerhet i lösningen av Navier-Stokes ekvationer av vätskeflöde, reduceras till över 16.000 kopplade och stela vanliga differentialekvationer. Slutligen diskuterar vi det breda spektrum av öppna forskningsteman som följer av det arbete som presenteras. Nyckelord: Bayesian Numerical Analysis, Osäkerhet Kvantifiering, Gaussian Processes, Differentialekvation Modeller, Osäkerhet i datormodeller. I alla vetenskaper, ekonomi och ingenjörskonst finns det en grundläggande tillit till användningen av differentialekvation modeller för att beskriva komplexa fenomen koncist, med hjälp av få men lätt att tolka parametrar, som vi betecknar med. I system av differentialekvationer, derivata med avseende på rumsliga variabler, x-D-D-D, och temporalvariabler, t-A-[a, b]-R+, är relaterade till implicit definierade tillstånd, u-x, t)-R-P, som därför ofta är analytiskt svåråtkomliga. Den största utmaningen med att arbeta med differentialekvationsmodeller, både ur matematiska och statistiska perspektiv, är att lösningarna i allmänhet inte är tillgängliga i sluten form. Följaktligen används numeriska lösare för att approximera en systemlösning över ett diskretiseringsgaller. Numerisk integration av ett system av differentialekvationer ger en deterministisk approximation baserad på subjektiva val avseende ordningen på den numeriska metoden, den angivna feltoleransen och det implicita diskretiseringsnätet, som alla påverkar kvaliteten på approximationen (Butcher, 2008; Higham, 1996). Dessa metoder gör implicit starka antaganden, till exempel att små förändringar i stegstorleken resulterar i små kumulativa förändringar i produktionen, och att det finns en unik, välutrustad lösning. Dessutom är meningsfulla slutsatser från modellen beroende av antagandet att de ungefärliga tillstånden innehåller försumbara numeriska fel. Även om den angivna feltoleransen och tillhörande punktvisa fel kan minskas genom att förfina diskretiseringsnätet, finns det naturligtvis en avräkning mellan den beräknade ansträngningen och noggrannheten i den approximerade lösningen, som avgör om antagandet av försumbart fel är rimligt i * oksana.chkrebtii@gmail.com; adress: Department of Statistics, The Ohio State University, 1958 Neil Avenue, 404 Cockins Hall, Columbus, OH 43210-1247 1 arXiv:1306365v2 [stat.ME] 28 Apr 2014 praxis. Lokala och globala numeriska fel definieras punktmässigt och relaterar till det asymptotiska beteendet hos modellens deterministiska approximation. Denna form av numerisk felanalys är inte lämplig för att kvantifiera den funktionella osäkerheten i lösningen i syfte att helt probabilistiska modellslutsatser. En ytterligare utmaning när man arbetar med modeller med differentialekvation är att klassiska approximationer av lösningar i vissa fall ger potentiellt vilseledande banor. Långtidsnumeriska lösningar kan vara globalt känsliga för trunkeringsfel som införs vid varje diskretiseringspunkt (t.ex. Sauer et al., 1997).............................................................. Illconditioned modeller kan ge upphov till en annan källa till osäkerhet i form av lösning multiplicity (Ascher et al., 1988; Beyn och Doedel, 1981). Detta kan ofta inte lätt verifieras analytiskt och utgör tydligt ett problem för klassiska numeriska integrationsmetoder som producerar en enda deterministisk lösning (Keller, 1968). I detta dokument tar vi upp de modelleringsutmaningar som uppstår vid användning av följande klasser av deterministiska vanliga och partiella differentialekvationsmodeller, för vilka analytiska lösningar i allmänhet inte finns tillgängliga. Vanliga differentialekvationsmodeller (ODE) representerar implicit derivatan av stater u(t), med avseende på tid t t genom ut (t) := dt u(t) = f (t, u(t), där vi gör standard antagandet att f är kontinuerlig i det första argumentet och Lipschitz kontinuerlig i det andra argumentet. Ingångar och gränsbegränsningar definierar flera modellvarianter. Initial Value Problem (IVP) modellerar systemet tillstånd med fast initialt tillstånd u * (a), utvecklas enligt ODE enligt följande, u t (t) = f t, u(t),  till, t [a, b], u(a) = u * (a). Förekomsten av en lösning garanteras under milda förhållanden (se till exempel Butcher, 2008; Coddington och Levinson, 1955). Sådana modeller kan vara högdimensionella och innehålla andra komplexiteter såsom algebraiska komponenter, funktionella ingångar eller högre ordertermer. Medan IVP-modeller anger ett fast initialt villkor på systemet tillstånd, kan Mixed Boundary Value Problem (MBVP) begränsa olika tillstånd vid olika tidpunkter. Vanligtvis är dessa begränsningar pålagda i slutet av tidsdomänen vilket ger den allmänna formen för två tillstånd blandade gränsvärdesproblem, som enkelt kan generaliseras till högre dimensioner och extrapoleras efter slutpunkten b. Medan en unik IVP-lösning existerar under relativt milda förhållanden, kan införandet av blandade gränsbegränsningar resultera i flera lösningar (t.ex. Keller, 1968) inför allvarliga problem för parameteruppskattningsmetoder. Den Delayy Initial Function Problem (DIFP) generaliserar den initiala begränsningen av IVPs till en initial funktion,  och därmed relatera derivatan av en process till både nuvarande och tidigare tillstånd vid släpar τ j  på [0, ∞), DIFPs är väl lämpade för att beskriva biologisk och fysisk dynamik som tar tid att föröka sig genom system. De utgör dock utmaningar för numeriska tekniker på grund av potentiellt stora och strukturerade trunkering fel (Bellen och Zennaro, 2003). Dessutom kan känsligheten hos DIFP-lösningar till små förändringar i den ursprungliga funktionen, såsom de på grund av interpoleringsfel, driva ett annars väluppfostrade system in i en kaotisk regim (Taylor och Campbell, 2007). Partiell differentialekvation (PDE) modeller representerar derivatan av tillstånd med avseende på flera argument, till exempel tid och rumsliga variabler. De huvudsakliga klasserna av PDE-modeller är baserade på elliptiska, paraboliska och hyperboliska ekvationer. Ytterligare öka deras komplexitet, funktionella gränsbegränsningar och inledande villkor gör PDE-modeller mer utmanande, medan deras underliggande teori är mindre utvecklad jämfört med ODE-modeller (Polyanin och Zaitsev, 2004 Numeriskt diskretiseringsfel karakteriseras ofta i form av en övre gräns, som kallas verifieringsfel i matematik och teknisk litteratur (Oberkampf och Roy, 2010). I de många fall där felet inte är försumbart, t.ex. i storskaliga Navier-Stokes-ekvationer som används vid modellering av vätskeflöde, krävs en korrekt representation av felutbredningen till följd av approximationen. Ett belysande exempel på kvantifiering av detta numeriska fel i fluid dynamics ges av Oliver et al. (2014) och ett preliminärt försök till en bayesisk formalisering av osäkerheten beskrivs av Oliver och Moser (2011). Det tillvägagångssätt som förespråkas i detta dokument är annorlunda till sin natur. Vi karakteriserar diskretiseringsosäkerhet med hjälp av en probabilistisk representation av lösningen beroende på fasta initiala förhållanden och lösningsparametrar. Vår metodik utmanar några av de antaganden som används av numeriska lösare och ger algoritmer robusta för många av de inferentiella problem som uppstår. Större uppmärksamhet ägnas åt de utmaningar som är förknippade med kvantifiering av osäkerheten i systemmodeller, särskilt de som bygger på matematiska beskrivningar av naturliga eller fysiska processer med hjälp av differentialekvationer (Ghanem och Spanos, 2003; Huttunen och Kaipio, 2007; Kaipio m.fl., 2004; Marzouk och Najm, 2009; Marzouk m.fl., 2007; Stuart, 2010). Den viktiga frågan om hur osäkerhet sprids genom komplexa system uppstår i en rad olika sammanhang, och fokus ligger ofta på två huvudproblem. Det framtida problemet handlar om hur indataosäkerheten (dvs. i parametrarna och de ursprungliga förhållandena) sprids över den numeriska banan. Denna typ av problem har undersökts i ingenjörsvetenskap och tillämpad matematik litteratur, som använder sig av provtagningsbaserade metoder, perturbering av initiala systemtillstånd (Mosbach och Turner, 2009), moment stängning metoder, och mer nyligen polynomiska kaos metoder (se Ghanem och Spanos, 2003; Xiu, 2009, för en fullständig översikt). Det omvända problemet gäller osäkerheten i de modellindata som getts mätningar med fel i modellens utgångar, och mycket arbete inom detta område har presenterats i statistiklitteraturen (Bock, 1983; Brunel, 2008; Calderhead och Girolami, 2011; Calderhead et al., 2009; Campbell and Lele, 2013; Campbell and Steele, 2011; Campbell and Chkrebtii, 2013; Dowd, 2007; Gugushvili och Klaasen, 2012; Ionides et al., 2006; Liang and Wu, 2008; Ramsay et al., 2007; Xue et al., 2010; Xun et al., 2014). Dessutom kan det omvända problemet i sig vara okonditionerat när data inte ligger i lösningsutrymmet för differentialekvationsmodellen (Brynjarsdóttir och O'Hagan, 2014; Kennedy och O'Hagan, 2001). Frågan om diskretiseringsosäkerhet är inneboende för både de framåtgående och omvänt problem, och det verkar naturligt att betrakta den okända lösningen av själva differentialekvationen som en del av den inferenta ramen. Det finns många exempel som tyder på att ett probabilistiskt funktionellt tillvägagångssätt är väl lämpat för att modellera lösningsosäkerhet för differentialekvationsmodeller; att lösa PDE:er innebär att uppskatta ett fält, DIFP:er kräver uppskattning av en okänd funktion, och MBVP:er kan ha flera lösningar. Alla dessa problem kan naturligt inramas i termer av åtgärder över ett funktionsutrymme (Stuart, 2010). Ett Bayesianskt funktionellt tillvägagångssätt erbjuder ett sätt att definiera och uppdatera åtgärder över möjliga banor för ett system med differentialekvationer och sprida den resulterande osäkerheten konsekvent genom den inferentiella processen. Den matematiska ramen för att definiera numeriska problem på funktionsytor har utvecklats genom det grundläggande arbetet Stuart (2010), Somersalo (2007), och Skilling (1991), och detta är den naturliga miljö som vi antar i vårt bidrag. Vi utvecklar dessa idéer ytterligare och undersöker allmänna tillvägagångssätt för att arbeta med klasser av differentialekvationer ur ett Bayesiskt perspektiv. Användningen av sannolikhet i numerisk analys har en lång historia. På 1950-talet utvecklades Monte Carlo metoder ursprungligen för att uppskatta analytiskt svåråtkomliga högdimensionella integraler, och förlängningar av denna metod har sedan dess haft en enorm inverkan på en mängd olika områden. Markovkedjan Monte Carlo (MCMC) metoder har gjort det möjligt att snabbt utveckla och tillämpa Bayesian statistiska metoder för att kvantifiera osäkerhet i vetenskapliga och tekniska problem. Teoretiska analyser av Bayesian statistiska metoder för klassisk numerisk analys problem verkar ha daterat hela vägen till Poincaré, som vältaligt sammanfattas av Diaconis (1988). På senare tid har O'Hagan (1992); Skilling (1991) också gett praktiska motiv för att utveckla och använda Bayesianska tekniker för att karakterisera den osäkerhet som är resultatet av beräkningsrestriktioner i fråga om integration, eftersom vi inte kan utvärdera en funktion vid varje punkt i dess inmatningsutrymme och därför måste ta hänsyn till diskretiseringsosäkerhet. Sådant arbete är nära relaterat till idén om dator emulering av stora och komplexa datormodeller; användningen av Gaussiska processer till modellfunktioner som är beräkning mycket dyra att utvärdera har funnit tillämpning inom många utmanande områden av vetenskap, såsom klimat förutsägelse och geofysik (Conti och O'Hagan, 2010; Kennedy 3 och O'Hagan, 2001; Tokmakian et al., 2012). Andra aktuella exempel på att tillämpa Bayesian strategier för numerisk integration inkluderar Hennig och Hauberg (2014), som tillämpar Skillings strategi för beräkning av Riemannian statistik, och Kennedy (1998); Osborne et al. (2012), vars metoder möjliggör en mer välinformerad probabilistisk strategi för valet av integrationspunkter, baserat på den förutsedda variationen i surrogatmodellen givet begränsade funktionsutvärderingar. Sådana metoder förutsätter dock inget fel i själva funktionsutvärderingen och är därför inte tillämpliga på modeller för differentialekvation, där funktionsutvärderingen endast ger ett ungefärligt derivat av systemets lösning vid nästa diskreta tidpunkt, baserat på de aktuella prediktiva fördelningarna över staterna. Syftet med detta dokument är att göra ett argument för att probabilistiskt karakterisera diskretiseringsfel vid lösning av allmänna klasser av differentialekvationer, och för att betrakta denna epistemiska osäkerhet som en viktig och integrerad del av den övergripande modellens inferensram. Vi utvecklar först en probabilistisk differentialekvationslösare som kännetecknar diskretiseringsosäkerhet i lösningen, genom att formalisera och avsevärt utvidga de idéer som först presenterades av Skilling (1991). Den sekventiella provtagningsmetod som föreslås är möjlig att beräkna genom rekursiv Bayesian uppdatering, och den följer noga den sekventiella konstruktionen av ett bevis på överensstämmelse som vi presenterar i tillägget, vilket säkerställer att algoritmen konvergerar till den exakta lösningen eftersom tidssteget tenderar att noll. I vårt tillvägagångssätt behandlar vi de sekventiella modellutvärderingarna som hjälpparametrar, över vilka lösningen sedan kan marginaliseras, och vi utvidgar den grundläggande metoden ytterligare till fall av problem med blandade gränsvärden som uppvisar multiplicitet av lösningar, försenar initiala funktionsproblem med delvis observerade initiala funktioner, stela och kaotiska PDEs via spektral projektion, samt ger en ram för att direkt lösa PDEs. Vi ger ett exempel på skalbarheten av vår strategi med hjälp av ett Navier-Stokes-system som omfattar lösningen på över 16.000 kopplade stela vanliga differentialekvationer. Slutligen antar vi en flexibel och allmän framåtsimuleringsmetod som gör det möjligt för oss att inbädda den formella kvantifieringen av diskretiseringsosäkerhet för differentialekvationer inom ramen för Kennedy och O'Hagan, och därmed definiera den som en del av det fullständiga inferenta förfarandet för inversa problem. Vi tillhandahåller kod för vår Probabilistic Differential Equation Solver (PODES), som tillåter replikering av alla resultat som presenteras i detta dokument. Detta finns på http://web.warwick.ac.uk/PODES. Vi börjar i avsnitt 2 med att diskutera modellavvikelser inom Bayesians strategi för kvantifiering av osäkerhet, som först presenterades av Kennedy och O'Hagan (2001) och därefter antogs i praktiken i hela statistiklitteraturen. Vi diskuterar hur befintliga numeriska metoder för integrering av differentialekvationer ignorerar diskretiseringsosäkerhet när lösningen inte är tillgänglig i sluten form. Vi utvecklar sedan en allmän Bayesiansk ram för att modellera denna epistemiska osäkerhet i differentialekvationer, inspirerad av arbetet av Skilling (1991) och O'Hagan (1992). Vi börjar i Avsnitt 3 med att definiera funktionella priors på stater och derivat för allmänna vanliga och partiella differentialekvationer. Sedan, i Avsnitt 4, utvecklar vi en Bayesian strategi för att uppdatera våra tidigare övertygelser givna modellinformation som erhållits sekventiellt och självkonsekvent över en given diskretisering rutnät, vilket probabilistiskt karakteriserar den finita dimensionella representationen av en underliggande infinit dimensional lösning givet modellparametrar och initiala villkor. Vi visar att den resulterande bakre banan med sannolikhet konvergerar till den exakta lösningen enligt standardantaganden. Vi införlivar sedan denna metod i det omvända problemet i avsnitt 5 och ger exempel i avsnitt 6. Slutligen diskuterar vi vårt föreslagna tillvägagångssätt för kvantifiering av Bayesisk osäkerhet för differentialekvationsmodeller och lyfter fram öppna forskningsområden. Tillägget innehåller bevis och ytterligare matematiska och algoritmiska detaljer. Vikten av att formellt kvantifiera källor till osäkerhet i datorsimuleringar när matematiskt modellera komplexa naturfenomen, såsom väder, havsströmmar, istäcke och cellulär proteintransport, är allmänt erkänd. Dessa modellbaserade simuleringar informerar resonemangsprocessen 4 när man t.ex. bedömer den ekonomiska risken vid beslut om oljefältshålskonfigurationer eller utformar en statlig politik som svar på extrema väderhändelser. Det är därför av stor betydelse att man tar hänsyn till alla osäkerhetskällor och sprider dem på ett enhetligt sätt under hela inferens- och beslutsprocessen. Överväg data, y(t), observeras vid diskreta tidpunkter t = [t 1, t 2,. .., t T ] och en uppsättning modellparametrar Med hjälp av den exakta lösningen av en matematisk modell som representeras av u * (t,  ska), en förenklad observationsmodell baserad på någon mätfelstruktur (t) är, I full allmänhet, en icke-linjär omvandling G u * (t, ) kan observeras, men för expositionell klarhet vi antar att detta är en identitet. I vår miljö, u * (t, ) är den unika funktionen som uppfyller ett allmänt system av differentialekvationer. När vi utför statistiska slutsatser över sådana modeller, är vi intresserade av den gemensamma bakre fördelningen över alla okända, eventuellt inklusive de initiala och gränsförhållanden. Genom hela papperet, kommer att utökas till att innehålla alla variabler av intresse. I Kennedys och O'Hagans banbrytande dokument (2001) identifierades ett antal osäkerhetskällor vid modellering av en naturlig eller fysisk process. Deras erkännande av osäkerhet på grund av ofullständig kunskap om processen och den motsvarande modellens otillräcklighet motiverar en probabilistisk syn på en missmatchning mellan modellerna δ(t) som också studerats i de senaste tidningarna av Brynjarsdóttir och O'Hagan (2014); Huttunen och Kaipio (2007); Kaipio och Solermo (2007); och Stuart (2010). I den tekniska litteraturen är denna form av fel också känd som valideringsfel (se t.ex. Oberkampf och Roy, 2010). Efter Kennedy och O'Hagan genom att definiera δ(t) som en slumpmässig funktion hämtad från en Gaussian Process (GP), blir observationsmodellen, På grund av bristen på en analytisk lösning för de flesta icke-linjära differentialekvationer, sannolikheten p y(t) u * (t,  på), kan inte uppnås i sluten form. Denna fråga behandlas i hela statistiklitteraturen genom att man ersätter den exakta sannolikheten med en surrogatlösning, p y(t)  på N (t, ) på grundval av en N-dimensionell ungefärlig lösning,û N (t, ) som erhållits med hjälp av numeriska integrationsmetoder (t.ex. en Runge-Kutta-lösare med N-tidsteg) vars noggrannhet har studerats väl (se t.ex. Henrici, 1964). Det finns dock fortfarande många scenarier för vilka standardmetoder för karakterisering av övre gränser för numeriska fel, och som därefter sprider denna osäkerhet under resten av den inferentiella processen är otillfredsställande. Begränsad beräkning och grov maskstorlek bidrar till det numeriska felet, som beskrivs i den omfattande översikten över Oberkampf och Roy (2010). Det till synes ofarliga antagandet av försumbara numeriska integrationsfel kan senare leda till allvarliga statistiska fördomar och vilseledande slutsatser för vissa klasser av differentialekvationsmodeller, och vi illustrerar detta i avsnitt 6.1 med ett enkelt exempel. Vi kan representera denna ytterligare osäkerhet genom att använda termen "T" = u Kennedy och O'Hagan (2001) föreslår modellering av datakoden med hjälp av en Gaussisk process som är agnostisk till den specifika formen av den underliggande matematiska modellen, som därför betraktas som en "svart låda". I diskussionen om deras papper hävdar H. Wynn att känslighetsekvationer och andra underliggande matematiska strukturer som styr datorsimuleringen också kan ingå i den övergripande osäkerhetsanalysen. I det aktuella dokumentet "öppnar vi den svarta lådan" genom att explicit modellera lösningen och tillhörande diskretiseringsosäkerhet,û N (t,  på) +. Detta gör att Kennedy och O'Hagan ramverket kan berikas ytterligare genom att införliva detaljerad kunskap om den matematiska modell som används. Observerar att den bakre mätningen över datormodellen definieras oberoende av de observerade data tyder på en ytterligare konstruktion genom att probabilistiskt karakterisera osäkerheten på grund av systemtillstånd som implicit definieras av differentialekvationer. Detta leder till utvecklingen av ett helt probabilistiskt system för att lösa differentialekvationer och föreslår en kraftfull ny metod för modellering av osäkerhet i datormodeller. Vi modellerar osäkerhet i en ändlig dimensionell representation av den oändliga dimensionslösningen genom en sannolikhetsangivelse på ett utrymme av lämpligt smidiga funktioner. Begränsa oss till Hilbert utrymmen för 5 modellera vår kunskap om u * (t,  på), definierar vi en Gaussian föregående åtgärd på funktionsutrymmet (Stuart, 2010). Vi modellerar sedan direkt vår kunskap om lösningen via den stokastiska processen u(t, ) och ersätter därmed (6) med y(t) = u(t,  ska) + δ(t) + (t). Bidraget i detta dokument är i vår definition och användning av u(t, ), därför för expositionell klarhet vi fokuserar vår uppmärksamhet på den gemensamma bakre mått över differentialekvation modell tillstånd, parametrar, och tillhörande hjälpparametrar,  på vår probabilistiska modell av osäkerhet, Vår Bayesian probabilistic integration ram för differentialekvation modeller utforskar fullt ut utrymmet för banor som ungefär uppfyller modelldynamik under teoretiska garantier för konsekvens. Den ram som presenteras i detta dokument demonstreras på modeller som beskrivs med hjälp av ickelinjära vanliga, fördröjningar, blandade gränsvärden och partiella differentialekvationer, och vi tar uttryckligen upp multiplicitet av lösningar i avsnitt 6.2. För notationell enkelhet vi nedan utelämnar beroendet av u på.... När man överväger osäkerhet i den oändliga dimensionella lösningen av differentialekvationer, är det naturliga måttet för utrymmen i en bred klass av funktioner Gaussian (Stuart, 2010). Som sådan definierar vi en Gaussian process (GP) föregående åtgärd gemensamt på staten u och dess tidsderivat u t. Vi kommer nu att diskutera den tidigare konstruktionen för vanliga differentialekvationer, och sedan överväga fallet för partiella differentialekvationer. I resten av tidningen kommer vi att beteckna med R λ en deterministisk, fyrkantig integrerbar kärna funktion med längdskala λ  och (0, ∞), och dess integrerade version av Q λ (t 1, t 2 ) = t1 a R λ (s, t 2 )ds. Här representerar a den nedre gränsen för vår tidsdomän. Vår modell för derivatet har covarians operator cov(, där α är en föregående precision parameter. Därför covariance på staten är dess integrerade version, cov(u(t 1 ), u(t 2 ) = α Korset covariance termer definieras på ett liknande sätt och betecknas som RQ(t 1, t 2 ) och QR(t 1, t 2 ) respektive. Vi noterar att RQ(t 1, t 2 ) = QR † (t 1, t 2 ), där † representerar angränsande. Vi antar en gemensam Gaussian föregående åtgärd på staten och dess derivat, Den Gaussiska villkorade tidigare åtgärder för t.ex.. p (u t (t)) u(t)) anta standardformulär (Stuart, 2010). I tillägget (avsnitt 8.2) diskuteras valet av tidigare medel och kovariansstruktur samt inverkan och valet av hjälpparametrar,  på = [α, λ]. För exposition, vi tar hjälpparametrarna  med fasta, även om vi noterar att propagera en distribution över  med hjälp av den probabilistiska lösaren är enkel, och uppskattningen av  och tas upp i avsnitt 6.4. Denna modell generaliserar enkelt till ODE problem av ordning större än en, genom att definiera den föregående gemensamt på stat, u, och alla nödvändiga derivat. Alternativt kan ODE-problem i högre ordning upprepas i den första ordningen genom att införa ytterligare stater. Vi har funnit att addition av tillstånd till modellen är beräkningsmässigt okomplicerat, så vi antar detta tillvägagångssätt i praktiken, arbetar med föregående i ekvation (8). PDE-problem kräver en föregående åtgärd över multivariata banor, samt modellering av enskilda eller högre orderderivat med avseende på rumsliga ingångar. Därför kommer den tidigare specifikationen att vara beroende av PDE-modellen. Betrakta som ett belysande exempel den paraboliska värmeekvationen, modellera värmespridningen över tiden längs en enda rumslig dimension av, 6 Ett modellering val för att införliva den rumsliga komponenten av PDE i covariansen är att anta en produktstruktur. Kovariansen över tiden definieras enligt ovan, och över rymden kan definieras med hjälp av en liknande konstruktion. Låt R μ vara en kärna funktion med längd skala μ. Q μ (z, x 2 )dz, där c betecknar en rumslig lägre gräns (i (9), detta är c = 0). De rumsliga kovariansstrukturerna definieras på liknande sätt som den timliga formen i ODE-exemplet som, där β är en rumslig föregående precisionsparameter. Kors covariances definieras analogt. Den tidigare konstruktionen följer genom att definiera en produktstruktur för rum och tid, I nästa avsnitt beskriver vi i detalj ramen för att sekventiellt koppla föregående i ekvation (8) med ODE-modellen. Vi beskriver motsvarande algoritm för uppdatering av föregående i ekvation (10) under värmeekvationen PDE i tillägget (avsnitt 8.1). Vi ger exempel som illustrerar dess användning experimentellt i Avsnitt 6. Den gemensamma föregående på staten och dess derivat definieras i ekvation (8) i termer av medelfunktioner m, m t och en kovariansstruktur som bestäms av valet av kärnfunktion R λ, som parameteriseras av hjälpvariablerna α och λ. Valet av kovarianskärnan bör återspegla våra antaganden om den exakta men okända differentialekvationslösningens smidighet. I tillägget (avsnitt 8.3) tillhandahåller vi covariansstrukturer baserade på två kärnor: den oändligt olika kvadratexponentiella och den icke-differentierbara enhetliga kärnan. Gaussiska processmodeller matchar vanligtvis kärnan för ett givet program till tidigare information om smidigheten i den underliggande funktionen (se t.ex. Rasmussen och Williams, 2006). Det är också viktigt att vara medveten om att införandet av orealistiskt strikta mjukhetsantaganden på statsrymden genom val av kovariansstruktur kan medföra en bedömningsfördel om den exakta lösningen inte är minst lika smidig. Därför, i fall där lösningen jämnhet inte är känd a priori, kan man fela på sidan av försiktighet genom att använda mindre vanliga kärnor. I de exempel som rapporteras i detta papper valde vi att arbeta med stationära derivata covariance strukturer för enkelhet, men vi påpekar att det kommer att finnas klasser av problem där icke-stationära kärnor kan vara lämpligare och kan också införlivas i den Gaussiska processramen. Medelfunktionen väljs utifrån alla förkunskaper om hur systemet utvecklas över tiden. Vanligtvis är dock sådan tidigare detaljerad kunskap inte tillgänglig, i vilket fall det är rimligt att använda en konstant medelfunktion. Dessutom, vi villkorar på de kända gränsvärdena, genomdriva dem genom att välja det föregående medelvärdet. Till exempel, för IVPs, vi uppfyller gränsen begränsningen exakt genom att välja en differentiable föregående genomsnittlig funktion m sådan att m(a) := u * (a). Nu när formen för den tidigare åtgärden har definierats, erhålls motsvarande bakre del i följande två avsnitt på ett sekventiellt sätt. Denna iterativa uppdatering följer noga den sekventiella strukturen för det bevis på överensstämmelse som anges i tillägget (avsnitt 8.6). Här påpekar vi att även om mängden och kvaliteten på tidigare information om den sanna lösningen kommer att påverka effektiviteten hos våra probabilistiska lösare, kan vi fortfarande förvänta oss konvergens eftersom tidssteget tenderar till noll, med förbehåll för vissa standardantaganden om kärnans funktion. Dessutom illustrerar vi i avsnitt 6.1 omfattningen av tidigare påverkan på den probabilistiska lösningen för värmeekvationen PDE. Vi anser att användningen av ett nollmedelvärde före på de rumsliga och temporal derivata har minimal inverkan på lösningen, även för någorlunda grov diskretisering rutnät. Vi inför nu en sekventiell ram för att karakterisera den epistemiska modellosäkerhetskomponenten i ekvation (7). Vi begränsar vår uppmärksamhet till ODE:s initialvärde (1) och överväger mer allmänna ODE- och PDE-problem i avsnitt 6. Vi modellerar diskretisering osäkerhet genom att uppdatera den gemensamma före på det okända tillståndet och derivatera iterativt, genom att utvärdera den deterministiska ODE-modellen vid ett ändligt antal N rutnät punkter över den temporal domänen, som vi betecknar med f 1:N. Algoritm 1 producerar ett prov från den gemensamma fördelningen, [ u, Låt oss överväga problemet med modellering osäkerhet i den exakta men okända lösningen av en icke-linjär ODE initialt värde problem. Vi ges differentialekvationsmodellen, u t (t) = f (t, u(t), som implicit definierar derivata i termer av en Lipschitz kontinuerlig funktion f av de stater som anges och ett fast initialvärde, u * (a). Vi modellerar vår osäkerhet om den unika explicita tillstånd, u, och dess derivat, u t, uppfyller det initiala värdeproblemet, via den gemensamma GP föregående åtgärd i ekvation (8). Följande förfarande länkar sekvensvis föregående på staten med ODE-modellen definieras av f. Överväga en diskretisering rutnät bestående av N tidspunkter Vi börjar med att fastställa det kända ursprungliga värdet, u(s 1 ) := u * (a) och beräkna den exakta derivatan f 1 := f (s 1, u(s 1 ),  på s 1, via den deterministiska ODE-modellen. Vi uppdaterar sedan vår gemensamma GP tidigare med tanke på den beräknade exakta derivata f 1, erhålla villkorlig prediktiv distribution för staten på den efterföljande rutnätet plats s 2, med, Denna prediktiva distribution beskriver vår nuvarande osäkerhet om lösningen vid tidpunkten s 2. Vi nu ta prov på ett förverkligande, u(s 2 ), av den prediktiva processen, och återigen länka vår före den deterministiska ODE-modellen genom att computing f 2 := f (s 2, u(s 2 ). I motsats till den första gången punkt s 1, Vi kan inte längre garantera att genomförandet vid den andra tidpunkten, u(s 2 ), och dess derivat, u t (s 2 ), exakt uppfyller ODE-modellen, dvs.. att u t (s 2 ) = f 2. Därför, vid tidpunkten s 2 vi modellerar uttryckligen missmatchningen mellan ODE utvärdering f 2 och processderivatet, u t (s 2 ), som, där omfattningen av missmatchningen kan beskrivas av variansen, C t (s 2, s 2 ), av den prediktiva bakre över derivatan ges av, För system där vi anser att denna missmatchning är starkt icke-gaussian, kan vi på lämpligt sätt ändra ovanstående modell och tillhörande provtagningsstrategi, som beskrivs i diskussionsavsnittet. Vi kan nu uppdatera vår gemensamma bakre från föregående iteration genom att konditionera på den utökade vektorn f 1:2 := [f 1, f (s 2, u(s 2 )]. Vi definierar därför matrisen ι 2×2 := diag{0, C t (s 2, s 2 )} för att beskriva missmatchningen mellan processderivatet och ODE-funktionen utvärderingar på s 1 och s 2. Den nya prediktiva bakre, Precis som i de föregående två stegen, vi prov på tillståndet förverkliga u (s 3 ) från ovanstående prediktiva bakre fördelning. Vi tillämpar sedan den deterministiska omvandlingen f, för att erhålla f 3 := f (s 3, u(s 3 ),  på vars missmatchning med det realiserade derivatet u t (s 3 ) åter modelleras som, Vi nästa förstärka f 1:3 := [f 1:2, f 3 ]. Motsvarande missmatchning matrisen har därför en diagonal struktur, ρ 3×3 := diag{ρ 2×2, C t (s 3, s 3 )}, där steg-förutsedd kovarians i derivatutrymmet är, De diagonala elementen i ρ, som är steg-förutsebara derivat varianser, är icke-decreasing, återspeglar vår växande osäkerhet om hur väl det realiserade tillståndet lyder ODE-modellen som vi tar vårt prov längre och längre bort från det kända ursprungliga tillståndet. Det visas också i tillägget (avsnitt 8.7) att dess beståndsdelar tenderar att nollställas med stegstorleken men i mycket snabbare takt än stegstorleken. Det allmänna systemet kan skrivas enligt algoritm 1, och en illustration av detta ges i figur 1. Algoritm 1 Prov från den gemensamma bakre distributionen av u och f 1:N för ett initialt ODE-värdesproblem givet..,., N Vid tidpunkten s 1 := a, initiera derivatet f 1 := f 1, u(s 1 ),... för initialtillstånd u(s 1 ) := u * (a) och definiera tillhörande modeller-derivativ missmatchning := 1×1 := 0; för n = 1 : N − 1 Definiera det prediktiva tillståndet medelvärde och varians, Prov steg-ahead realisering u(s n+1 ) från den prediktiva fördelningen av tillståndet, utvärdera ODE-modellen f n+1 := f (s n+1, u(s n+1 ),...) för realisering u(s n+1 ) vid efterföljande rutnät punkt, s n+1, och förstärka vektorn f 1:n+1 := [f 1:n +1 ]; Definiera den prediktiva derivatvariansen, och öka matrisen Vi betonar att integreringen av ODE-modellen fortsätter probiliserande via Gaus integration. Vi får en bakre fördelning över banor 9 Figur 1 : Illustration av Algoritm 1 för att generera ett prov från den gemensamma fördelningen av derivat observationer och möjliga banor med densitet p u(t), f 1:N. Med tanke på två derivatmodellers realiserande (röda punkter) får vi en bakre fördelning över derivatutrymmet (överst till vänster) och över tillståndsutrymmet (nedre mitten). Ett prov dras sedan från den prediktiva bakre delen över staterna vid nästa tidpunkt s n (nedre mitten), och en modell realiseras genom att kartlägga u (s n ) till derivatutrymmet via funktionen f (topp mitt, högra röd punkt). Med tanke på dessa tre modellutvärderingar kan detta förfarande upprepas (nederst till höger, överst till höger) på ett liknande sätt. 10 u som styrs av differentialekvationsmodellen över ett diskret antal rutnätspunkter. I tillägget (avsnitt 8.6) ger vi ett bevis på att den bakre process u(t) som erhålls via Algoritm 1 under vissa förhållanden tenderar att vara den unika lösning u * (t) som uppfyller IVP (1), enligt följande, där h är den maximala längden mellan på varandra följande diskretiseringspunkter. Vi presenterar nu en ram för att sprida modellens diskretiseringsosäkerhet som karakteriseras i avsnitt 4 genom det omvända problemet, där vi vill dra slutsatser om modellparametrar från mätdata. Vi beskriver en möjlig Markov kedja Monte Carlo förfarande som genererar ett prov från den gemensamma bakre fördelningen av tillståndet på dataplatserna t = [t 1, · · · ·, t T ] och de okända parametrarna  på villkor att bullriga observationer y(t) = G u(t), + (t) av staterna. Vi antar att parametrarna  och inkluderar alla okända initiala villkor, eller hjälpparametrar. För enkelhetens skull fokuserar vi återigen på direkta observationer av stater som styrs av ett initialt värdeproblem. Utvidgningen till andra ODE- och PDE-problem är dock okomplicerad med en probabilistisk lösning. Vi tillhandahåller en applikation där stater indirekt observeras genom den icke-linjära omvandlingen av stater G i avsnitt 6.4. Algoritm 2 mål den bakre densiteten i ekvation (7) genom framåt modellförslag via villkorlig simulering från Algoritm 1. Genom detta steg i förslaget undviker man behovet av att uttryckligen beräkna den indragbara marginaltätheten p u(t), f 1:N  på ett sådant sätt att den kan genomföras på ett effektivt sätt som beskrivs i bilagan (avsnitt 8.1). Sådana delvis sannolikhetsfria MCMC-implementationer (Marjoram et al., 2003) används i stor utsträckning inom området inferens för stokastiska differentialekvationer (se exempelvis Golightly och Wilkinson, 2011) för att simulera provvägar inom det omvända problemet. Algoritm 2 Rita K-prover från den bakre fördelningen med densitet p  för, u(t) y(t),  på ett villkorligt prov en realiserande av tillståndet u(t) via Algoritm 1; för k = 1 : K gör Föreslå  en  på förslagstäthet, där q är en förslagstäthet; Prov ett probabilistiskt förverkligande av tillståndet u (t) betingat på  ska via Algoritm 1; Compute: Vi har nu alla komponenter att ta en helt Bayesian metod för att kvantifiera osäkerhet på differentialekvation modeller av naturliga och fysiska system. I det här avsnittet använder vi det ramverk som utvecklats i de två föregående avsnitten i tillämpningar till ett brett spektrum av system, inklusive problem med ODE- och PDE-gränsvärde och fördröjda initiala funktionsproblem. 11 Som ett belysande exempel visar vi användningen av vår probabilistiska ram för värmeekvationen som presenteras i (9). Vi modellerar vår osäkerhet om lösningen genom den tidigare (10) definierade över tid och rum med hjälp av en produkts kovariansstruktur som beskrivs i avsnitt 3. Vi använder vårt probabilistiska ramverk för att integrera tillståndet över tiden vid var och en av de rumsliga diskretiseringsgridpunkterna, med användning av en enhetlig kovarianskärna för tidskomponenten och en kvadratexponentiell kovarianskärna för den rumsliga komponenten. Vi karakteriserar därför diskretiseringsosäkerheten i både tid och rumsliga domäner genom konditionering av PDE-modellutvärderingar som motsvarar ungefärliga tidsderivat vid varje diskretiseringsgridplats, beroende på sekventiellt urval av realiserande från den prediktiva bakre fördelningen över den andra ordningens rumsliga derivat. Vi beskriver den fullständiga algoritmiska konstruktionen i tillägget (Algoritm 3). I följande numeriska simuleringar tar vi hänsyn till dynamiken med σ = 1 och initial funktion u * (x, 0) = sin (xπ), x  till [0, 1]. Vi betraktar först det probabilistiska framtidsproblemet med hjälp av en mängd olika diskretiseringsgaller. I figur 2 jämför vi den exakta lösningen med den probabilistiska lösningen på två olika rutnät; en grov diskretisering av 15 punkter i den rumsliga domänen och 50 punkter i den tidsmässiga domänen, och en finare diskretisering av 29 punkter i den rumsliga domänen och 100 punkter i den tidsmässiga domänen. Vi observerar från simuleringarna att i takt med att maskstorleken blir finare minskar osäkerheten i lösningen, i överensstämmelse med konsekvensresultatet (bilaga 8.6), där det visas att den probabilistiska lösningen bör tendera till den exakta lösningen eftersom gridavståndet tenderar att vara noll. Även om detta är ett belysande exempel med hjälp av ett enkelt leksakssystem, är karakteriseringen av rumslig osäkerhet avgörande för mer komplexa modeller, där det finns beräkningsbegränsningar som begränsar antalet systemutvärderingar som kan utföras. Vi kan se effekten av sådan diskretiseringsosäkerhet genom att utföra bakre inferens för parametern σ given data simulerad från en exakt lösning med σ = 1. Figur 3 visar den bakre fördelningen jämfört med den som erhålls genom att använda både en "framåt i tiden, centrerad i rymden" (FTCS) ändlig skillnadslösare och en probabilistisk lösare under en mängd olika diskretiseringsgaller. Användningen av en deterministisk lösare illustrerar problemet med inferential bias och överkonfident bakre varians som kan uppstå om diskretisering osäkerhet inte beaktas och alltför grov ett rutnät används. I denna illustrativa miljö förökar användningen av en probabilistisk lösning osäkerheten om diskretisering i lösningen fram till den bakre fördelningen över parametrarna. Vi får parameteruppskattningar som tilldelar positiv sannolikhetsmassa till det sanna värdet av δ, även när man använder ett grovt diskretiserat rutnät, undvika problemet med överkonfident parameterslutsatser som utesluter det verkliga värdet. MBVP introducerar utmaningar för många befintliga numeriska lösare, som vanligtvis förlitar sig på optimering och teorin om IVP för att uppskatta det ospecificerade initialtillståndet, u(a) i (2). Optimering över u(a) utförs tills motsvarande IVP-lösning uppfyller det angivna gränstillståndet u * (b) inom en användarspecifik tolerans. För expositionell enkelhet, anser vi ett gränsvärde problem med en begränsning belägen vid varje gräns för domänen, nämligen v(a), u(b) = v * (a), u * (b). Vi behandlar gränsvärdet u * (b) som en datapunkt och betraktar lösningen som ett inferensproblem över det ospecificerade initialvärdet u(a). Sannolikheten definierar därför missmatchningen mellan gränsvärdet u(b), som erhålls från den realiserade probabilistiska lösningen med en del u(a), med det exakta gränsvärdet u * (b), enligt följande, där m (u) (b) och C (u) (b) b) är det bakre medelvärdet och kovariansen för tillstånd u vid tidpunkten b erhålls från utvärderingen av Algoritm 1. Den bakre fördelningen av staterna har därför densitet, som uppvisar multimodalitet över staterna, som visas på vänster sida av figur 4. Medan deterministiska numeriska lösare förlitar sig på en tolerans för ad hoc-slutpunkters missmatchning, är den probabilistiska ramen naturligt figur 2: Vi illustrerar den probabilistiska produktionen av lösningen till värmeekvationen PDE, med σ = 1, integrerad mellan t = 0 och t = 0,25 med hjälp av två rutnätsstorlekar; den grövre maskan (visas i blått) består av 15 rumsliga diskretiseringspunkter och 50 tids diskretiseringspunkter, den finare maskan består av 29 rumsliga diskretiseringspunkter och 100 tids diskretiseringspunkter. Vi visar de rumsliga bakre förutsägelserna vid tre tidpunkter; t = 0,02 (överst), t = 0,12 (mellan) och t = 0,22 (underst). Den exakta lösningen vid varje tidpunkt representeras av den gröna linjen. Felstängerna visar medelvärden och 2 standardavvikelser för var och en av de probabilistiska lösningarna beräknade med hjälp av 50 simuleringar. Figur 3: Vi illustrerar det omvända problemet genom att dra slutsatser om parametern σ i värmeekvationen, integrerad mellan t = 0 och t = 0,25. Vi genererar data över ett rutnät med 8 rumsliga diskretiseringspunkter och 25 diskretiseringspunkter genom att använda den exakta lösningen med σ = 1, och sedan addera buller med standardavvikelsen 0,005. Vi använder först den probabilistiska differentialekvationslösaren (PODES) med tre rutnätsstorlekar; en grov maska bestående av 8 rumsliga diskretiseringspunkter och 25 diskretiseringspunkter (långt till vänster), en finare maska bestående av 15 rumsliga diskretiseringspunkter och 50 diskretiseringspunkter (andra från vänster) och ytterligare en finare maska bestående av 29 rumsliga diskretiseringspunkter och 100 diskretiseringspunkter (andra från höger). Observera förändringen i skala när den bakre variansen minskar med ökande upplösning av diskretiseringen. Som en belysande jämförelse visar vi de bakre distributionerna med hjälp av en deterministisk framåt i tiden, centrerad i rymden (FTCS) integrationssystem (långt till höger). Om diskretiseringen inte är tillräckligt bra, får vi en överkonfidentiell fördomsfull bakdel som tilldelar försumbar sannolikhetsmassa till det verkliga värdet av δ. Däremot ger användningen av den exakta lösningen, som förväntat, en fullständigt objektiv baksida. Figur 4: 3000 prov togs från den bakre sannolikhetstätheten (13) och banorna visas för båda tillstånden (vänster, över och under). Den marginella bakre sannolikhetstätheten över det okända initiala tillståndet u( 1 2 ) visas också (höger). definierar denna tolerans genom en prediktiv fördelning av (12). En viktig faktor för att lösa MBVPs probabilistically är att prov effektivt från potentiellt multimodala posteriors över det initiala värdet. Vi rekommenderar därför användning av en lämplig MCMC-system som parallell härdning (Geyer, 1991), som snabbt kan identifiera och utforska olika regioner med hög bakre sannolikhet. Vi tillhandahåller ett sådant genomförande av parallell härdning i tillägget (Algoritm 5). Som en demonstration av probabilistiskt lösa ett blandat gränsvärde problem, anser vi ett specialfall Lane-Emden modellen, som används för att beskriva densitet u av gasformiga sfäriska objekt, såsom stjärnor, som en funktion av dess radie t, (Shampine, 2003). Vi skriver om den kanoniska andra ordningen ODE som ett system av en första ordningen ekvationer med ett blandat gränsvärde, med gränsförhållanden, u * (1) = ε 3/2 och v * ( 1 2 ) = −288/2197. Det okända initialtillståndet u * ( 1 2 ) tilldelas en diffus Gaussian tidigare med medelvärde 1,5 och standardavvikelse 2 u * (1) - v * ( 1 2 ), vilket återspeglar möjligheten att flera lösningar kan finnas i ett brett spektrum av initialtillstånd. I detta exempel valde vi den fyrkantiga exponentiella kovariansen för att modellera vad vi förväntar oss att vara en mycket smidig lösning. Diskretiseringsnätet består av 100 jämnt fördelade punkter. Längdskalan är inställd på två gånger stegstorleken för diskretiseringsgaller och den tidigare precisionen är inställd på 1. Figur 4 visar ett bakre prov från ekvation 13 och identifierar två regioner med hög densitet i det bakre som motsvarar distinkta banor som ungefär uppfyller modelldynamiken med tanke på vårt diskretiseringsgaller. Figur 4 visar också multimodalitet genom den marginella bakre över det okända initialtillståndet u( 1 2 ). Detta exempel illustrerar behovet av att korrekt modellera lösningsosäkerheten på ett funktionellt sätt. En numerisk approximation, även med motsvarande numeriska felgränser, kommer inte att upptäcka förekomsten av en andra lösning, medan en probabilistisk strategi gör det möjligt för oss att bestämma antal och plats för möjliga lösningar. I detta fall är en probabilistisk metod för att lösa problemet med differentialekvation inriktad på modellfördomar som orsakas av många olika lösningar. 14 Vi presenterar följande exempel som ett bevis på konceptet att den probabilistiska lösaren kan tillämpas tillförlitligt och enkelt på ett mycket högdimensionellt dynamiskt system. Systemet Navier-Stokes är en grundläggande modell av flytande dynamik, som införlivar lagar om bevarande av massa, energi och linjärt momentum, samt fysiska egenskaper hos en okomprimerbar vätska över vissa domäner givna begränsningar som påtvingas längs gränserna. Det är en viktig komponent i komplexa modeller inom oceanografi, väder, luftföroreningar och glaciärrörelser. Trots dess omfattande användning, dynamiken i Navier-Stokes modeller är dåligt förstådd även på små tidsskalor, där de kan ge upphov till turbulens. Vi anser Navier-Stokes PDE-modellen för tidsutveckling av 2 komponenter av hastigheten, u : D → R 2, av en okomprimerbar vätska på en torus, D := [0, 2π, ) × [0, 2π], uttryckt i sfäriska koordinater. Gränsvärdesproblemet Navier-Stokes definieras av: där  på: = är den Laplacian operatören sådan att u = u x2x2, och är lutningsoperatören sådan att u = u x1 + u (15) inte är känd i sluten form. Ofta är mängden av intresse den lokala spinning rörelse av den okomprimerbara vätska, som kallas vorticitet, som vi definierar som, = − till × u, där  till × u representerar rotationscurl definieras som korsprodukt av  på och u, med positiv vorticitet motsvarande medurs rotation. Denna variabel kommer att användas för att bättre visualisera den probabilistiska lösningen av Navier-Stokes systemet genom att reducera de två komponenterna av hastigheten till en dimensionell funktion. Vi diskretiserar Navier-Stokes-modellen (15) över ett rutnät av storlek 128 i varje rumslig dimension. Därför, en pseudo spektral projektion i Fourier utrymme ger 16,384 kopplade, stela ODS med tillhörande begränsningar. Fullständiga uppgifter om pseudospektral projektion ges för att möjliggöra fullständig replikering av dessa resultat och finns på den åtföljande webbplatsen. Figur 6.3 visar fyra framåt simulerade vorticitetsbanor (långa rader), som erhålls från två hastighetskomponenter som styrs av Navier-Stokes ekvationer (15) vid fyra olika tidpunkter (långa kolumner). Lätta skillnader i tillståndsdynamiken kan ses vid den sista tidpunkten, där de fyra banorna tydligt avviker från varandra. Dessa skillnader uttrycker den epistemiska osäkerhet som uppstår till följd av diskretisering av den exakta men okända oändligt dimensionella lösningen. Det probabilistiska tillvägagångssättet för framåtsimulering av DIFP beskrivs i Algoritm 6 i tillägget. Förutom att modellera diskretisering osäkerhet på ett strukturerat, funktionellt sätt, vårt ramverk gör det möjligt för oss att enkelt införliva osäkerhet i den ursprungliga funktionen genom framåt simulering; faktiskt, den ursprungliga funktionen kan det sig själv endast vara tillgänglig vid ett begränsat antal noder. Det probabilistiska tillvägagångssättet kvantifierar den osäkerhet som är förenad med uppskattningen av (t) och sprider den rekursivt genom staterna. Även när  och t är helt specificerade i förväg påverkar det nuvarande tillståndets beroende av tidigare beräknade tillstånd noggrannheten hos den numeriska lösningen med hjälp av standardlösare, även på kort sikt. I följande exempel redogör vi för osäkerheten i samband med aktuella och försenade uppskattningar av systemtillstånd, och tar helt itu med dessa potentiella biaskällor i uppskattningen av modellparametrar. 16 deras ursprungliga tillstånd, återvänder till cytoplasman för att användas i nästa aktiveringscykel. Detta sista steg är inte väl förstådd och är proxied i modellen av den okända tidsfördröjningen τ. Modellen för denna mekanism beskriver förändringar i 4 reaktionstillstånd hos STAT-5 genom den ickelinjära DIFP, De initiala funktionskomponenterna  och (2) (t) =  och (3) (t) =. (4) (t) är överallt noll, medan den konstanta initialfunktionen  och (1) (t) är okänd. Tillstånden för detta system kan inte mätas direkt, men observeras genom en icke-linjär omvandling, och parameteriseras av de okända skalningsfaktorerna k 5 och k 6. De indirekta mätningarna av staterna antas vara förorenade med additivt noll-medelvärde Gaussian buller, ε(t), med experimentellt fastställda standardavvikelser, Vår analys är baserad på experimentella data som uppmätts på platserna, t, från Swameye et al. (2003), som består av 16 mätningar för de två första tillstånden i observationsprocessen. Raue m.fl. (2009) Utnyttja ytterligare en artificiell datapunkt för var och en av de tredje och fjärde observationsprocesserna för att hantera bristen på parameteridentifierbarhet för detta system. Därför antar vi detta antagande i vår analys. Drivfunktionen, EpoR A : [0, 60] → R +, är inte känd, men mätt vid 16 diskreta tidpunkter t EP O. Enligt Raue et al. (2009), Vi antar att observationer mäts utan fel. Vi utgår vidare från att denna funktion har samma smidighet som lösningstillståndet (delvis linjärt förstaderivat). Den fullständiga villkorliga fördelningen av forceringfunktionen ges av en GP interpolation av observationer. Forward inference för denna modell kommer att användas inom det statistiska inversa problemet med att återvinna okända parametrar och första starttillstånd,  på = [k 1,. .., k 6, τ, u (1) (0)] från experimentella data. Vi demonstrerar fullt probabilistic inference för stat banor och parametrar av den utmanande 4 tillstånd fördröjning initial funktionsmodell (16) som beskriver dynamiken i JAK-STAT cellulär signal transduction bana (Raue et al., 2009 ). Det har gjorts flera analyser av Jak-STAT väg mekanism baserat på dessa data (t.ex. Campbell och Chkrebtii, 2013; Raue et al., 2009; Schmidl et al., 2003; Swameye et al., 2003). Trots de olika modelleringsantaganden som olika författare har beaktat, liksom olika slutledningsmetoder, har några intressanta gemensamma drag identifierats som motiverar den explicita modelleringen av diskretiseringsosäkerhet för denna ansökan. För det första har de inkorrekta och beräkningsmässiga begränsningar av numeriska tekniker som krävs för att lösa systemekvationerna fått vissa författare att ta till grova ODE- approximationer eller till och med indirekt kringgå numeriska lösningar via Generalised Smoothing. Detta motiverar en formell analys av diskretiseringsfelets struktur och utbredning genom det omvända problemet. En annan fråga är att modellen (16) och dess varianter är missriktade. Ovanstående studier 17 tyder på att modellen inte är tillräckligt flexibel för att passa de tillgängliga uppgifterna, men det är inte klart hur mycket av detta fel som beror på modellavvikelser och hur mycket som kan tillskrivas diskretiseringsfel. Vår analys fortskrider genom att definiera tidigare fördelningar på de okända parametrarna enligt följande, Vi fick prover från den bakre fördelningen av modellparametrarna,  och = [k 1,. .., k 6, τ, u (1) (0)], lösningstillstånd, u(t,  på), och hjälpvariabler,  på grund av data y(t) med MCMC. För att konstruera en Markov-kedja som effektivt färdas genom parameterutrymmet för denna multimodala bakre fördelning använde vi en parallell härdningsprovtagare (Geyer, 1991) med 10 parallella kedjor längs en jämnt fördelade temperaturprofil över intervallet [0,5, 1]. Varje probabilistisk DIFP-simulering genererades med hjälp av algoritm 6 under ett lika stort rutnät för diskretisering av storlek N = 500. Fullständiga algoritmiska detaljer finns i tillägget (Algoritm 8). Vi fick två grupper av bakre prover, var och en av storlek 50 000. Inom kedjan bedömdes konvergensen genom att man testade lika metoder mellan olika iterationsintervaller i kedjan (Geweke 1992). Mellan kedjekonvergensen bedömdes på liknande sätt. Dessutom såg vi till att acceptfrekvensen föll ungefär inom det accepterade intervallet på 18-28 % för var och en av de två parameterblocken, och att den totala acceptfrekvensen för förflyttningar mellan två kedjor låg kvar på ungefär 5–15 %. Korrelationsytor och beräkningar av kärnans densitet med föregående för marginalparameterns bakre del visas i figur 6. Alla parametrar utom den tidigare precisionen α identifieras med hjälp av data, inklusive hastighetsparametern k 2 som endast tycks vara svagt identifierad. Vi observerar starka korrelationer mellan parametrar, i enlighet med tidigare studier på detta system. Till exempel finns det en stark korrelation mellan skalningsparametrarna k 5, k 6 och det första tillståndet u (1) (0). Intressant nog verkar det finnas ett samband mellan den probabilistiska lösarens längdskala λ och den första, tredje och fjärde reaktionsfrekvensen. Dessutom har denna korrelation en icke-linjär struktur, där den högsta parameterdensitetsregionen tycks förändras med längdskala som innebär stark känslighet för lösningsspecifikationerna. Längdskalan är den probabilistiska analogen, under en begränsad kovarians, av stegnumret i en numerisk metod. Men i analyser baserade på numerisk integration, valet av en numerisk teknik effektivt fixerar denna parameter till ett godtyckligt värde som väljs a priori. Vårt resultat här tyder på att för detta problem är de infererade parametervärdena starkt och icke-linjärt beroende av valet av den numeriska metod som används. Vi spekulerar om att denna effekt kan bli ganska allvarlig för mer känsliga system, så att ignorera diskretisering osäkerhet inom det omvända problemet kan leda till inferentiella fördomar. Ett prov från den marginella bakre delen av statsbanorna och motsvarande observationsprocess visas i figur 7. Felraderna på datapunkterna visar två standardavvikelser för mätfelet från Swameye et al. (2003).................................................................. Det står omedelbart klart att vår modell, som innehåller diskretiseringsosäkerhet i det framtida problemet, fortfarande inte helt fångar dynamiken i de observerade uppgifterna. Denna systematiska brist på passform tyder på att det finns skillnader i modellen utöver den som beskrivs av diskretiseringsosäkerheten. Detta dokument har presenterat en probabilistisk formalism för att beskriva strukturen av ungefärlig lösningsosäkerhet för allmänna system av differentialekvationer. I stället för att tillhandahålla en enda uppsättning diskreta funktionsvärden som ungefär motsvarar de begränsningar som systemet med differentialekvationer medför ger vårt tillvägagångssätt en sannolikhetsmått över utrymmet för sådana oändliga dimensionella funktioner. Detta är en avvikelse : Marginal bakre fördelning av modellparametrarna med hjälp av probabilistisk framåt simulering baserad på ett prov av storlek 50 000, genereras med hjälp av en parallell härdning algoritm med tio kedjor. Sannolikheten för täthet visas i svart. 19 från den befintliga vedertagna praxisen att använda en deterministisk numerisk integrationskod för att erhålla en ungefärlig statistisk felmodell som en del av slutsatssprocessen. Detta gör det möjligt att anta sannolikhetskalkyler för att på ett enhetligt sätt sprida funktionell osäkerhet när man löser det statistiska inversa problemet. Vi har visat att vårt probabilistiska ramverk redan kan tillämpas på komplexa och högdimensionella system med hjälp av exemplet Navier-Stokes. Det krävs dock ytterligare arbete för att förbättra beräkningseffektiviteten till en punkt där den är jämförbar med befintliga genomföranden av numeriska integrationskoder. I princip, eftersom den probabilistiska integrationsmetoden vi har föreslagit bygger på lösningar av linjära system, delar den samma algoritmiska skalning som de flesta implicita numeriska integrationskoder (t.ex. Crank-Nicholson), och man räknar med att algoritmisk och kodutveckling med tiden kommer att säkerställa att vår metodik uppnår liknande prestandanivåer och blir en del av standardverktygslådan för klimatforskare, geoforskare och ingenjörer. Algoritmiska framsteg och kodutveckling kommer att vara ett spännande och givande område för pågående utredningar med stor genomslagskraft. Utökning av metoden till system som uppvisar riktningsderivat fel är omedelbar, till exempel i fluid dynamics applikationer. I avsnitt 4 har vi faktiskt modellerat missmatchningen mellan föregående på derivatan och ODE-modellen som Gaussian. Avslappning av detta antagande garanterar inte längre en sluten form representation av den uppdaterade tidigare på statsderivatet, men detta kan övervinnas genom ytterligare ett skikt av Monte Carlo provtagning inom Algoritm 1. Sådana frågor skulle dessutom motivera utvecklingen av effektiva MCMC-algoritmer i detta sammanhang. I detta arbete har vi direkt kunnat utnyttja den matematiska modellens struktur inom inferensprocessen på ett mer informerat sätt än genom att behandla det oåtkomliga systemet med differentialekvationer och tillhörande numerisk lösningskod som en "svart låda". Vi kan också utnyttja den inbyggda mångfaldiga struktur som induceras av den probabilistiska modellen senare för att öka effektiviteten i en mängd olika MCMC-metoder (Girolami och Calderhead, 2011). Som vi redan har noterat kan kontrollen av de hjälpparametrar som är förknippade med jämnhet i utrymmet för möjliga banor också ge flexibilitet till MCMC provtagningsmetoder. Detta inslag i vår probabilistiska metod gör det möjligt att definiera mellanliggande måltätheter för population MCMC-metoder såsom Smooth Functional Tempering (Campbell and Steele, 2011), vilket ökar effektiviteten när man utforskar de komplexa bakre densiteter som uppstår i tillämpningar som beskrivs av system med ickelinjära vanliga eller partiella differentialekvationer. Den föreslagna probabilistiska metoden tar formen av en linjär funktionell projektion, vilket gör det möjligt att direkt beräkna känsligheten i systemtillstånden med avseende på parametrar. Lokal känslighetsanalys kan vara användbar i tekniska tillämpningar, slutsatser, eller för att underlätta optimering. Detta skulle också göra det möjligt att till exempel införliva vår probabilistiska ram i en klassisk icke-linjär ram med minst kvadrater (Bates and Watts, 2007). Det probabilistiska tillvägagångssätt som vi förespråkar gör det i grund och botten möjligt att definiera en formell avvägning mellan osäkerhet som orsakas av numerisk noggrannhet och det omvända problemets storlek, genom val av diskretiseringsnät. Detta problem är av stort intresse för osäkerhet kvantifiering gemenskapen (se t.ex. Arridge et al., 2006; Kaipio et al., 2004). Vi kan nu kvantifiera och jämföra de enskilda bidragen till den övergripande osäkerheten från diskretisering, modellfelspecifikation och Monte Carlo-fel inom det probabilistiska ramverket. Vår probabilistiska modell av diskretiseringsosäkerhet kan användas för att styra meshförfining och till och med meshdesign för komplexa modeller, eftersom en prediktiv fördelning nu finns tillgänglig, vilket utgör grunden för experimentutformningsmetoder för meshutveckling. Numeriska lösningar av ODE och PDE för komplexa modeller är också beroende av adaptivt val av mesh. Den probabilistiska formalism som presenteras här kan nu informera mesh val på ett probabilistiskt sätt, genom att optimera ett valt konstruktionskriterium. Vissa arbeten i denna riktning har redan utförts i Chkrebtii (2013), där det naturliga Kullback-Leibler divergenskriteriet med framgång används för att adaptivt välja diskretiseringsnätet probabilistiskt. Ett sådant tillvägagångssätt föreslogs också av Skilling (1991) med hjälp av korsentropin. Slutsats och förutsägelse för datorexperiment (Sacks et al., 1989 ) bygger på numeriska lösningar av storskaliga systemmodeller. För närvarande är numerisk osäkerhet i modellen till stor del ignoreras, även om det informellt införlivas genom en covariance nugget i emulatorn (se, till exempel, Gramacy och Lee, 2012). Att anta en probabilistisk strategi i stor skala kommer att få praktiska konsekvenser på detta område genom att göra det möjligt att lätta på det felfria antagande som antas när man modellerar utdata från datorkoder, vilket leder till mer realistiska och flexibla emulatorer. Vid modellering av osäkerhet om en exakt men okänd lösning kan den probabilistiska integrationsformalismen anses ge både en uppskattning av den exakta lösningen och dess tillhörande funktionsfelanalys. Många befintliga numeriska metoder kan faktiskt tolkas inom ramen för denna allmänna probabilistiska ram. Detta tyder på möjligheten att både generalisera befintliga numeriska lösare och att utveckla nya provtagningssystem för att komplettera de probabilistiska lösningar som vi beskriver i detta bidrag. Kaotiska system uppstår i modellering av ett stort antal fysiska system, inklusive laserhålor, kemiska reaktioner, vätskerörelse, kristalltillväxt, väderprognoser och jordbävningsdynamik (Baker och Gollub, 1996, kapitel 7). Extrem känslighet för små störningar kännetecknar kaotisk dynamik, där effekten av diskretiseringsosäkerhet blir ett viktigt bidrag till global lösningsosäkerhet. Även om systemets långsiktiga beteende är förankrat i de ursprungliga tillstånden, resulterar förekomsten av numeriska diskretiseringsfel i exponentiella avvikelser från den exakta lösningen. Följaktligen sönderfaller information om initialtillstånd snabbt i takt med att systemlösningen utvecklas i tid. Denna insikt, demonstreras och förklaras av (Berliner, 1991), visas i följande exempel. Vi anser den klassiska Lorenz initialvärde problem (Lorenz, 1963), en bedrägligt enkel tre-state ODE modell av konvektiv vätskerörelse framkallad av en temperaturskillnad mellan en övre och nedre yta. Ett prov från detta systems probabilistiska bakre del visas i figur 8 givna modellparametrar i det kaotiska systemet och ett fast initialtillstånd. Det finns en kort tidsrymd inom vilken det finns försumbar osäkerhet i lösningen, men ackumuleringen av osäkerheten resulterar snabbt i divergerande, men ändå mycket strukturerade, lösningar eftersom flödet begränsas kring en kaotisk attraherare. Detta exempel belyser behovet av en funktionell modell av diskretiseringsosäkerhet för att ersätta numeriska punktvisa felgränser. Inom det omvända problemet gör ett probabilistiskt tillvägagångssätt också att vi kan kvantifiera förlusten av modellinformation när vi går bort från det ursprungliga tillståndet. Studien av kaotiska system är helt klart ett mycket viktigt forskningsområde, och man räknar med att det probabilistiska tillvägagångssätt som föreslås kan bidra till att påskynda framstegen. Vi återkommer nu till en viktig diskussion om den klass av problem för vilka probabilistisk integration är väl lämpad. Som vi har sett med slutledning för JAK-STAT-modellen, även relativt liten diskretisering osäkerhet kan förstärkas genom den icke-linjära framåt-modellen och kan införa skevhet i uppskattningen av parametrar som banan är mycket känslig för. Man kan dock anta att slutsatser om stabila system med låga dimensioner, utan topologiska begränsningar av lösningen, inte får ha några betydande diskretiseringseffekter. Fördelarna med vår formulering ligger dock utan tvekan vid forskningens spets när det gäller kvantifiering av osäkerhet, att hantera massiva ickelinjära modeller, att uppvisa komplex eller till och med kaotisk dynamik, och starka rumsliga-geometriska effekter (t.ex. modeller för undervattensflöde). Att lösa sådana system på ett ungefär är ett problem som fortfarande är i utkanten av aktuell forskning inom numerisk analys. Vi föreslår att ett probabilistiskt tillvägagångssätt skulle ge en viktig ytterligare uppsättning verktyg för att studera sådana system. Sammanfattningsvis hoppas vi att vår uppsats kommer att uppmuntra mycket mer forskning vid detta spännande gränssnitt mellan matematisk analys och statistisk metodik.
Ett liknande tillvägagångssätt tillämpas i REF, vilket ger en icke-parametrisk bakre snarare än en Gaussian approximation.
Bayesian Solution Uncertainty Quantification for Differential Equations
{'venue': None, 'journal': 'arXiv: Methodology', 'mag_field_of_study': ['Mathematics']}