690 likes | 888 Views
Intr o duc t i on. Computatio n in fields. Accelerating fields. Conclusion s and F uture W o rk. References. Ap p endix. Accelerating the ‘ fields ’ P ac k age in R : The o ry and Application John P aige 1 ( paigejo@gmail.com )
E N D
Introduction Computationinfields Acceleratingfields ConclusionsandFutureWork References Appendix Acceleratingthe‘fields’PackageinR: TheoryandApplication JohnPaige1(paigejo@gmail.com) IsaacLyngaas2(irlyngaas@gmail.com) SrinathVadlamani3(srinathv@ucar.edu) DougNychka3(nychka@ucar.edu) 1MacalesterCollege 2FloridaStateUniversity 3NationalCenterofAtmosphericResearch July31,2014 tudentswhoare
Introduction Computationinfields Acceleratingfields ConclusionsandFutureWork References Appendix Introduction IntroductiontoProblemGoalsandMotivation Computationinfields MathematicalFoundation Acceleratingfields EigenandCholeskyDecompositionsParameterOptimization ConclusionsandFutureWork Appendix
Introduction Computationinfields Acceleratingfields ConclusionsandFutureWork References Appendix AnalyzingSpatialDatawithfields:Kriging Colorado Average Spring High Temperature (Celsius) 41 20 40 15 Latitude 39 10 38 5 37 −108 −106 −104 −102 Longitude Datasetfromfieldspackage(Nychkaetal.,2014b)
Introduction Computationinfields Acceleratingfields ConclusionsandFutureWork References Appendix AnalyzingSpatialDatawithfields:Kriging Colorado Average Spring High Temperature (Celsius) 41 20 40 ? 15 Latitude 39 10 38 5 37 −108 −106 −104 −102 Longitude Datasetfromfieldspackage(Nychkaetal.,2014b)
Introduction Computationinfields Acceleratingfields ConclusionsandFutureWork References Appendix AnalyzingSpatialDatawithfields:Kriging Datasetfromfieldspackage(Nychkaetal.,2014b)Krigingestimateforλ=.103,θ=5.45
Introduction Computationinfields Acceleratingfields ConclusionsandFutureWork References Appendix DifficulttoAnalyzeLargeDatasets • DifficultiesinKriging: • Krigingwithmanyobservations Fields' CO2Dataset 50 Latitude 0 −50 −150 −100−50050 100 150 Longitude CO2datasetinfieldspackage(Nychkaetal.,2014b)
Introduction Computationinfields Acceleratingfields ConclusionsandFutureWork References Appendix DifficulttoAnalyzeLargeDatasets • DifficultiesinKriging: • Krigingwithmanyobservations Fields' CO2Dataset 50 Latitude 0 −50 −150 −100−50050 100 150 Longitude CO2datasetinfieldspackage(Nychkaetal.,2014b)
Introduction Computationinfields Acceleratingfields ConclusionsandFutureWork References Appendix DifficulttoAnalyzeLargeDatasets • DifficultiesinKriging: • Krigingwithmanyobservations • • O(n3) Fields' CO2Dataset 50 Latitude 0 −50 −150 −100−50050 100 150 Longitude CO2datasetinfieldspackage(Nychkaetal.,2014b)
Introduction Computationinfields Acceleratingfields ConclusionsandFutureWork References Appendix DifficulttoAnalyzeLargeDatasets • DifficultiesinKriging: • Krigingwithmanyobservations • • O(n3) • 100,000observations,fixedparameters:∼8hours Fields' CO2Dataset 50 Latitude 0 −50 −150 −100−50050 100 150 Longitude CO2datasetinfieldspackage(Nychkaetal.,2014b)
Introduction Computationinfields Acceleratingfields ConclusionsandFutureWork References Appendix DifficulttoAnalyzeLargeDatasets • DifficultiesinKriging: • Krigingwithmanyobservations • • O(n3) • 100,000observations,fixedparameters:∼8hours • 100,000observations,10 parametersamples:∼80hours Fields' CO2Dataset 50 Latitude 0 −50 −150 −100−50050100150 Longitude CO2datasetinfieldspackage(Nychkaetal.,2014b)
Introduction Computationinfields Acceleratingfields ConclusionsandFutureWork References Appendix DifficulttoAnalyzeLargeDatasets • DifficultiesinKriging: • Krigingwithmanyobservations • • O(n3) • 100,000observations,fixedparameters:∼8hours • 100,000observations,10 parametersamples:∼80hours • Krigingovertimeformanyparametersets: Fields' CO2Dataset 50 Latitude 0 −50 −150 −100−50050100150 Longitude CO2datasetinfieldspackage(Nychkaetal.,2014b)
Introduction Computationinfields Acceleratingfields ConclusionsandFutureWork References Appendix DifficulttoAnalyzeLargeDatasets • DifficultiesinKriging: • Krigingwithmanyobservations • • O(n3) • 100,000observations,fixedparameters:∼8hours • 100,000observations,10 parametersamples:∼80hours • Krigingovertimeformanyparametersets: • 5,000observations,30years(monthlydata),20parametersamples:nearly14hours Fields' CO2Dataset 50 Latitude 0 −50 −150 −100−50050100150 Longitude CO2datasetinfieldspackage(Nychkaetal.,2014b)
Introduction Computationinfields Acceleratingfields ConclusionsandFutureWork References Appendix ResearchGoals • LearnKrigingtheory • AcceleratethefieldspackageinR • Choleskyandeigendecompositions • Spatialparameterestimation(maximumlikelihoodestimation)
Introduction Computationinfields Acceleratingfields ConclusionsandFutureWork References Appendix Why Focusonfields? • Spatial Problem Workflow Times • goeR • fields(withmKrig) • MadeinR • (Nychkaetal.,2014b) • Free,opensource • Popularwithstatisticians • Easytouse 100 Time(Minutes) 10 1 0.1 0 2000 4000 6000 NumberofObservations 8000 10000 geoRandfieldspackagesusedinthisplot:DiggleandRibeiro(2007);Nychkaetal.(2014b);RibeiroJrandDiggle(2001)
Introduction Computationinfields Acceleratingfields ConclusionsandFutureWork References Appendix Why Focusonfields? • Spatial Problem Workflow Times • goeR • fields(withmKrig) • MadeinR • (Nychkaetal.,2014b) • Free,opensource • Popularwithstatisticians • Easytouse • Fast(forR) 100 Time(Minutes) 10 1 0.1 0 2000 4000 6000 NumberofObservations 8000 10000 geoRandfieldspackagesusedinthisplot:DiggleandRibeiro(2007);Nychkaetal.(2014b);RibeiroJrandDiggle(2001)
Introduction Computationinfields Acceleratingfields ConclusionsandFutureWork References Appendix • Why FocusonMatrixDecompositionsandOptimization? • EigenandCholeskydecompositionstakealongtime • Over10,000observationsmeansover40%ofthecomputationtime Spatial Surface Estimation Computation Time Complete Workflow Workflow CholeskyDecompositions 0 1 2 3 4 5 6 Time(Minutes) • Spatialparameteroptimization requiresmatrixdecompositions • Betteroptimizationmeans fewermatrixdecompositions 0 2000 4000 6000 Number of Observations 8000 10000 • Note: runningonaCalderanodeon Jellystone • GPUs:2TeslaM2070-Q • CPUs:28-core2.6-GhzIntelXeonE5-2670(SandyBridge)
Introduction Computationinfields Acceleratingfields ConclusionsandFutureWork References Appendix • Why FocusonMatrixDecompositionsandOptimization? • EigenandCholeskydecompositionstakealongtime • Over10,000observationsmeansover40%ofthecomputationtime Percent Workflow Time Taken by Cholesky Decompositions 20 40 60 80 100 • Spatialparameteroptimization requiresmatrixdecompositions • Betteroptimizationmeans fewermatrixdecompositions Percent 0 0 2000 4000 6000 Number of Observations 8000 10000 • Note: runningonaCalderanodeon Jellystone • GPUs:2TeslaM2070-Q • CPUs:28-core2.6-GhzIntelXeonE5-2670(SandyBridge)
Introduction Computationinfields Acceleratingfields ConclusionsandFutureWork References Appendix • Why FocusonMatrixDecompositionsandOptimization? • EigenandCholeskydecompositionstakealongtime • Over10,000observationsmeansover40%ofthecomputationtime Percent Workflow Time Taken by Cholesky Decompositions 20 40 60 80 100 • Spatialparameteroptimization requiresmatrixdecompositions • Betteroptimizationmeans fewermatrixdecompositions Percent 0 0 2000 4000 6000 Number of Observations 8000 10000 • Note: runningonaCalderanodeon Jellystone • GPUs:2TeslaM2070-Q • CPUs:28-core2.6-GhzIntelXeonE5-2670(SandyBridge)
Introduction Computationinfields Acceleratingfields ConclusionsandFutureWork References Appendix • Why FocusonMatrixDecompositionsandOptimization? • EigenandCholeskydecompositionstakealongtime • Over10,000observationsmeansover40%ofthecomputationtime Percent Workflow Time Taken by Cholesky Decompositions 20 40 60 80 100 • Spatialparameteroptimization requiresmatrixdecompositions • Betteroptimizationmeans fewermatrixdecompositions Percent 0 0 2000 4000 6000 Number of Observations 8000 10000 • Note: runningonaCalderanodeon Jellystone • GPUs:2TeslaM2070-Q • CPUs:28-core2.6-GhzIntelXeonE5-2670(SandyBridge)
Introduction Computationinfields Acceleratingfields ConclusionsandFutureWork References Appendix • Why FocusonMatrixDecompositionsandOptimization? • EigenandCholeskydecompositionstakealongtime • Over10,000observationsmeansover40%ofthecomputationtime Percent Workflow Time Taken by Cholesky Decompositions 20 40 60 80 100 • Spatialparameteroptimization requiresmatrixdecompositions • Betteroptimizationmeans fewermatrixdecompositions Percent 0 0 2000 4000 6000 Number of Observations 8000 10000 • Note: runningonaCalderanodeon Jellystone • GPUs:2TeslaM2070-Q • CPUs:28-core2.6-GhzIntelXeonE5-2670(SandyBridge)
IntroductionComputationinfields Acceleratingfields ConclusionsandFutureWork References Appendix KrigingSet-Up(WithNoTrendComponent) yi =g(xi)+ei yi:observations xi:locations g:Krigingsurface ei:errorfollowingN(0,σ2) distribution,independent
IntroductionComputationinfields Acceleratingfields ConclusionsandFutureWork References Appendix KrigingSet-Up(WithNoTrendComponent) yi =g(xi)+ei yi:observations xi:locations g:Krigingsurface ei:errorfollowingN(0,σ2) distribution,independent
IntroductionComputationinfields Acceleratingfields ConclusionsandFutureWork References Appendix KrigingSet-Up(WithNoTrendComponent) yi =g(xi)+ei yi:observations xi:locations g:Krigingsurface ei:errorfollowingN(0,σ2) distribution,independent
IntroductionComputationinfields Acceleratingfields ConclusionsandFutureWork References Appendix KrigingSet-Up(WithNoTrendComponent) yi =g(xi)+ei yi:observations xi:locations g:Krigingsurface ei:errorfollowingN(0,σ2) distribution,independent
IntroductionComputationinfields Acceleratingfields ConclusionsandFutureWork References Appendix KrigingSet-Up(WithNoTrendComponent) yi =g(xi)+ei yi:observations xi:locations g:Krigingsurface ei:errorfollowingN(0,σ2) distribution,independent
IntroductionComputationinfields Acceleratingfields ConclusionsandFutureWork References Appendix • KrigingSet-Up(WithNoTrendComponent) • yi =g(xi)+ei • yi:observations • xi:locations • g:Krigingsurface • ei:errorfollowingN(0,σ2) distribution,independentAssumptions(Nychkaetal.,2014a): • • E(g(x))=0 • Krigingsurfacehaszeromean • • E[g(x)g(xt)]=ρ·k(x,xt): • k:correlationfunction (e.g.e−||x−xt||/θ) • ρ:correlationstrength(orsignalstrength)
IntroductionComputationinfields Acceleratingfields ConclusionsandFutureWork References Appendix • KrigingSet-Up(WithNoTrendComponent) • yi =g(xi)+ei • yi:observations • xi:locations • g:Krigingsurface • ei:errorfollowingN(0,σ2) distribution,independentAssumptions(Nychkaetal.,2014a): • • E(g(x))=0 • Krigingsurfacehaszeromean • • E[g(x)g(xt)]=ρ·k(x,xt): • k:correlationfunction (e.g.e−||x−xt||/θ) • ρ:correlationstrength(orsignalstrength)
IntroductionComputationinfields Acceleratingfields ConclusionsandFutureWork References Appendix • KrigingSet-Up(WithNoTrendComponent) • yi =g(xi)+ei • yi:observations • xi:locations • g:Krigingsurface • ei:errorfollowingN(0,σ2) distribution,independentAssumptions(Nychkaetal.,2014a): • • E(g(x))=0 • Krigingsurfacehaszeromean • • E[g(x)g(xt)]=ρ·k(x,xt): • k:correlationfunction (e.g.e−||x−xt||/θ) • ρ:correlationstrength(orsignalstrength)
IntroductionComputationinfields Acceleratingfields ConclusionsandFutureWork References Appendix • KrigingSet-Up(WithNoTrendComponent) • yi =g(xi)+ei • yi:observations • xi:locations • g:Krigingsurface • ei:errorfollowingN(0,σ2) distribution,independentAssumptions(Nychkaetal.,2014a): • • E(g(x))=0 • Krigingsurfacehaszeromean • • E[g(x)g(xt)]=ρ·k(x,xt): • k:correlationfunction (e.g.e−||x−xt||/θ) • ρ:correlationstrength(orsignalstrength) • Reparameterization: λ=σ2/ρ:smoothingparameter (1/λ: signal/noiseratio)
Introduction Computationinfields Acceleratingfields ConclusionsandFutureWork References Appendix ParameterOptimization:MaximumLikelihoodEstimation Definition Alikelihoodfunction,L,givesthechancethatasetofobservationsoccurinamodelgiventhemodelparameters. Definition αˆMLEisamaximumlikelihoodestimateofαifαˆMLEmaximizes thedatalikelihood: L(αˆMLE)≥L(α),∀α
Introduction Computationinfields Acceleratingfields ConclusionsandFutureWork References Appendix ParameterOptimization:MaximumLikelihoodEstimation Definition Alikelihoodfunction,L,givesthechancethatasetofobservationsoccurinamodelgiventhemodelparameters. Definition αˆMLEisamaximumlikelihoodestimateofαifαˆMLEmaximizes thedatalikelihood: L(αˆMLE)≥L(α),∀α
Introduction Computationinfields Acceleratingfields ConclusionsandFutureWork References Appendix ParameterOptimization:MaximumLikelihoodEstimation Definition Alikelihoodfunction,L,givesthechancethatasetofobservationsoccurinamodelgiventhemodelparameters. Definition αˆMLEisamaximumlikelihoodestimateofαifαˆMLEmaximizes thedatalikelihood: L(αˆMLE)≥L(α),∀α
Introduction Computationinfields Acceleratingfields ConclusionsandFutureWork References Appendix ParameterOptimization:MaximumLikelihoodEstimation Definition Alikelihoodfunction,L,givesthechancethatasetofobservationsoccurinamodelgiventhemodelparameters. Definition αˆMLEisamaximumlikelihoodestimateofαifαˆMLEmaximizes thedatalikelihood: L(αˆMLE)≥L(α),∀α Thedatalog-likelihoodgiveninNychkaetal.(2014a)foragiven θandλis: yT(Σθ+λI)−1y1 lnL(θ,λ)=− + ln|Σθ+λI|+C 2 2
Introduction Computationinfields Acceleratingfields ConclusionsandFutureWork References Appendix ParameterOptimizationinfields yT(Σθ+λI)−1y1 lnL(θ,λ)=− + ln|Σθ+λI|+C 2 2 • y:vectorofobservationvalues • Σθ:covariancematrix,where: • (Σθ)ij=ρk(xi,xj)=ρe−||xi−xj||/θ • C:constant
Introduction Computationinfields Acceleratingfields ConclusionsandFutureWork References Appendix ParameterOptimizationinfields yT(Σθ+λI)−1y1 lnL(θ,λ)=− + ln|Σθ+λI|+C 2 2 • y:vectorofobservationvalues • Σθ:covariancematrix,where: • (Σθ)ij=ρk(xi,xj)=ρe−||xi−xj||/θ • C:constant
Introduction Computationinfields Acceleratingfields ConclusionsandFutureWork References Appendix ParameterOptimizationinfields yT(Σθ+λI)−1y1 lnL(θ,λ)=− + ln|Σθ+λI|+C 2 2 • y:vectorofobservationvalues • Σθ:covariancematrix,where: • (Σθ)ij=ρk(xi,xj)=ρe−||xi−xj||/θ • C:constant
Introduction Computationinfields Acceleratingfields ConclusionsandFutureWork References Appendix ParameterOptimizationinfields yT(Σθ+λI)−1y1 lnL(θ,λ)=− + ln|Σθ+λI|+C 2 2 • y:vectorofobservationvalues • Σθ:covariancematrix,where: • (Σθ)ij=ρk(xi,xj)=ρe−||xi−xj||/θ • C:constant
Introduction Computationinfields Acceleratingfields ConclusionsandFutureWork References Appendix ParameterOptimizationinfields yT(Σθ+λI)−1y1 lnL(θ,λ)=− + ln|Σθ+λI|+C 2 2 • y:vectorofobservationvalues • Σθ:covariancematrix,where: • (Σθ)ij=ρk(xi,xj)=ρe−||xi−xj||/θ • C:constant • Howcanthislikelihoodbecomputedquickly?
Introduction Computationinfields Acceleratingfields ConclusionsandFutureWork References Appendix ParameterOptimizationinfields yT(Σθ+λI)−1y1 lnL(θ,λ)=− + ln|Σθ+λI|+C 2 2 • y:vectorofobservationvalues • Σθ:covariancematrix,where: • (Σθ)ij=ρk(xi,xj)=ρe−||xi−xj||/θ • C:constant • Howcanthislikelihoodbecomputedquickly? • Kriguseseigendecomposition • mKrigusesCholeskydecomposition
Introduction Computationinfields Acceleratingfields ConclusionsandFutureWork References Appendix ParameterOptimizationinfields yT(Σθ+λI)−1y1 lnL(θ,λ)=− + ln|Σθ+λI|+C 2 2 • y:vectorofobservationvalues • Σθ:covariancematrix,where: • (Σθ)ij=ρk(xi,xj)=ρe−||xi−xj||/θ • C:constant • Howcanthislikelihoodbecomputedquickly? • Kriguseseigendecomposition • mKrigusesCholeskydecomposition
Introduction Computationinfields Acceleratingfields ConclusionsandFutureWork References Appendix Acceleratingfields:EigenandCholeskyDecompositions • MAGMA(Agulloetal.,2009) • Freelyavailablelibrarywithmulti-GPUcomputingcapability • MuchfasterthandefaultR
Introduction Computationinfields Acceleratingfields ConclusionsandFutureWork References Appendix AcceleratedWorkflowTime Spatial Problem Workflow Times Default R One GPU Two GPUs 8 Time (Minutes) 6 4 2 0 0 5000 10000 15000 Number of Observations • >2500observations⇒ acceleratedworkflowisfaster • >10000observations⇒ ≥1.55×speedup(for1or2GPUs)
Introduction Computationinfields Acceleratingfields ConclusionsandFutureWork References Appendix AcceleratedWorkflowTime Spatial Problem Workflow Times Default R One GPU Two GPUs 8 Time (Minutes) 6 4 2 0 0 5000 10000 15000 Number of Observations • >13000observations⇒ • twoGPUworkflowfasterthanoneGPUby≥1sec
Introduction Computationinfields Acceleratingfields ConclusionsandFutureWork References Appendix Acceleratingfields:MaximumLikelihoodEstimation • Goal: • Quicklyfindsetofmodel parametersmaximizingdata likelihood • Questions: • IsEigenorCholeskydecompositionfasterformaximizinglikelihood? • Howdodifferentwaysofsplitting upthesedecompositionsamongcoresandGPUsaffectlikelihood maximizationtime? ObservationsfromCO2datasetinfieldspackage(Nychkaetal., 2014b)
Introduction Computationinfields Acceleratingfields ConclusionsandFutureWork References Appendix Acceleratingfields:MaximumLikelihoodEstimation • Goal: • Quicklyfindsetofmodel parametersmaximizingdata likelihood • Questions: • IsEigenorCholeskydecompositionfasterformaximizinglikelihood? • Howdodifferentwaysofsplitting upthesedecompositionsamongcoresandGPUsaffectlikelihood maximizationtime? ObservationsfromCO2datasetinfieldspackage(Nychkaetal., 2014b)
Introduction Computationinfields Acceleratingfields ConclusionsandFutureWork References Appendix Acceleratingfields:MaximumLikelihoodEstimation • Goal: • Quicklyfindsetofmodel parametersmaximizingdata likelihood • Questions: • IsEigenorCholeskydecompositionfasterformaximizinglikelihood? • Howdodifferentwaysofsplitting upthesedecompositionsamongcoresandGPUsaffectlikelihood maximizationtime? ObservationsfromCO2datasetinfieldspackage(Nychkaetal., 2014b)
Introduction Computationinfields Acceleratingfields ConclusionsandFutureWork References Appendix Acceleratingfields:MaximumLikelihoodEstimation • Goal: • Quicklyfindsetofmodel parametersmaximizingdata likelihood • Questions: • IsEigenorCholeskydecompositionfasterformaximizinglikelihood? • Howdodifferentwaysofsplitting upthesedecompositionsamongcoresandGPUsaffectlikelihood maximizationtime? ObservationsfromCO2datasetinfieldspackage(Nychkaetal., 2014b)
Introduction Computationinfields Acceleratingfields ConclusionsandFutureWork References Appendix CholeskyWinsOverEigendecomposition Cholesky Decomposition Speedup Over Eigendecomposition (Default Implementations) 0 15 40 60 80 100 120 140 Speedup 0 5000 10000 15000 NumberofObservations • Optimizingλforfixedθ: • 10to15Choleskydecompositions • Oneeigendecomposition
Introduction Computationinfields Acceleratingfields ConclusionsandFutureWork References Appendix CholeskyWinsOverEigendecomposition Cholesky Decomposition Speedup Over Eigendecomposition (Default Implementations) 0 15 40 60 80 100 120 140 Speedup 0 5000 10000 15000 NumberofObservations • <20000observations⇒ Choleskyisatleast18×faster
Introduction Computationinfields Acceleratingfields ConclusionsandFutureWork References Appendix CholeskyWinsOverEigendecomposition Cholesky Decomposition Speedup Over Eigendecomposition (Default Implementations) 0 15 40 60 80 100 120 140 Speedup 0 5000 10000 15000 NumberofObservations • CholeskydecompositionachievedbetterspeedupswithGPUs