Learn R Programming

GPoM (version 1.0)

gloMoId: gloMoId : global Model Identification

Description

Algorithm for global modelling in a Polynomial and canonical formulation of Ordinary Differential Equations. Univariate Global modeling aims to obtain multidimensional models from single time series (Gouesbet & Letellier 1994, Mangiarotti et al. 2012). An example of such application can be found in Mangiarotti et al. (2014) For a multivariate application, see GPoMo (Mangiarotti 2015, Mangiarotti et al. 2016).

Example: For a model dimension nVar=3, the global model will read: \(dX1/dt = X2\) \(dX2/dt = X3\) \(dX3/dt = P(X1,X2,X3).\)

Usage

gloMoId(series, tin = NULL, dt = NULL, nVar = NULL, dMax = 1,
  weight = NULL, show = 1, filterReg = NULL, winL = 9)

Arguments

series

The original data set: either a single vector corresponding to the original variable; Or a matrix containing the original variable in the first column and its successive derivatives in the next columns. In the latter case, for the construction of n-dimensional model, series should have \(nVar+1\) columns since one more derivative will be necessary to identify the model parameters. Variable nVar will be set equal to n. In the former case, that is when only a single vector is provided, the derivatives will be automatically recomputed. Therefore, the dimension nVar expected for the model has to be provided.

tin

Input date vector which length should correspond to the variables of the input data (same number of lines).

dt

The time sampling of the input series.

nVar

The model dimension expected. This parameter will be deduced from the input data (series) if series is a matrix. If series is a vector, the expected dimension nVar should be provided.

dMax

Maximum degree of the polynomial functions allowed in the model (see poLabs).

weight

Weighting function of input data series. By default uniform weight is applied. This weighting function can also be used to apply the analysis piecewise using zeros and ones.

show

Indicates (2) or not (0-1) the algorithm progress

filterReg

A vector that specifies the a priori structure of the the function of the output model in canonical form. The organisation of the filter follows the convention defined in poLabs and can be obtained using poLabs(nVar,dMax). Value is 1 if the regressor is available, 0 if it is not.

winL

Total number of points used for computing the derivatives of the input time series (data). This parameter will be used as an input in function drvSucc.

Value

A list of five elements :

init Initial time series and the successive derivatives used in the modeling.

filterReg Structure of the output model. Value is 1 if the regressor is available, 0 if it is not. The corresponding terms can be obtained with poLabs(dMax,dMax).

codereconstr Simulated trajectory: the initial variable and its derivatives obtained by integrating the equations of the identified model.

codeK Values of the identified coefficients corresponding to the regressors defined in filterReg.

coderesTot Residual signal of the model.

coderesSsMod Residual signal of the closer submodels.

References

[1] Gouesbet G., Letellier C., Global vector-field reconstruction by using a multivariate polynomial L2 approximation on nets, Physical Review E, 49 (6), 4955-4972, 1994. [2] Mangiarotti S., Coudret R., Drapeau L., & Jarlan L., Polynomial search and global modeling : Two algorithms for modeling chaos, Physical Review E, 86, 046205, 2012. [3] Mangiarotti S., Drapeau L. & Letellier C., Two chaotic models for cereal crops observed from satellite in northern Morocco. Chaos, 24(2), 023130, 2014. [4] Mangiarotti S., Low dimensional chaotic models for the plague epidemic in Bombay (18961911), Chaos, Solitons & Fractals, 81(A), 184-196, 2015. [5] Mangiarotti S., Peyre M. & Huc M., A chaotic model for the epidemic of Ebola Virus Disease in West Africa (2013-2016). Chaos, 26, 113112, 2016.

See Also

gPoMo, autoGPoMoSearch, autoGPoMoTest, poLabs

Examples

Run this code
# NOT RUN {
#############
# Example 1 #
#############
# load data
#Example 2
data("Ross76")
tin <- Ross76[,1]
data <- Ross76[,2:3]

# Application with a simplified structure
reg <- gloMoId(data[0:500,2], dt=1/100, nVar=2, dMax=2, show=0)

#############
# Example 2 #
#############
# load data
data(NDVI)

# Definition of the Model structure
terms <- c(1, 0, 0, 0, 1, 1, 1, 1, 0, 1, 1, 1, 0, 1, 1, 1, 1, 1, 1, 1)
poLabs(3,3,terms==1)
reg <- gloMoId(NDVI [,1:1], dt=1/125, nVar=3, dMax=3,
               show=0, filterReg=terms==1)


# }
# NOT RUN {
#############
# Example 3 #
#############
# load data
data("Ross76")
# time vector
tin <- Ross76[1:500,1]
# single time series
series <- Ross76[1:500,3]
# some noise is added
series[1:100] <- series[1:100] + 0.01 * runif(1:100, min = -1, max = 1)
series[301:320] <- series[301:320] + 0.05 * runif(1:20, min = -1, max = 1)
# weighting function
W <- tin * 0 + 1
W[1:100] <- 0  # the first fourty values will not be considered
W[301:320] <- 0  # twenty other values will not be considered either
reg <- gloMoId(series, dt=1/100, weight = W, nVar=3, dMax=2, show=1)
visuEq(3, 2, reg$K, approx = 4)
# first weight which value not equal to zero:
i1 = which(reg$finalWeight == 1)[1]
v0 <-  reg$init[i1,1:3]

reconstr <- numicano(nVar=3, dMax=2, Istep=5000, onestep=1/250, PolyTerms=reg$K,
                     v0=v0, method="ode45")
plot(reconstr$reconstr[,2], reconstr$reconstr[,3], type='l', lwd = 3,
                            main='phase portrait', xlab='time t', ylab = 'x(t)', col='orange')
# original data:
lines(reg$init[,1], reg$init[,2], type='l',
      main='phase portrait', xlab='x', ylab = 'dx/dt', col='black')
# initial condition
lines(v0[1], v0[2], type = 'p', col = 'red')

# }

Run the code above in your browser using DataLab