ajustement des paramètres des EDO lors de l'utilisation du solveur d'EDO octave/matlab

Question

J'utilise OdePkg dans Octave pour résoudre un système d'EDO rigides, par exemple par ode5r :

function yprime = myODEs(t,Y,param)
    yprime = [
        - param(1) * Y(1);                      # ODE for Y(1)
        param(1) * Y(1) - param(2) Y(2) * Y(3); # ODE for Y(2)
        param(2) Y(2) * Y(3)                    # ODE for Y(3)
                                                # etc.
];

time_span = [1, 24]         # time span of interest
Y0        = [1.0, 1.1, 1.3] # initial values for dependent variables Y
param     = [7.2, 8.6, 9.5] # parameters, to be optimized

[t, Y] = ode5r(@myODEs, time_span, Y0, ..., param);

Le solveur stocke les variables dépendantes Y dans une matrice par rapport au temps t (vecteur) :

t     Y(1)  Y(2)  Y(3)
0.0   1.0   1.1   1.3
0.1   ...   ...   ...
0.5   ...   ...   ...
0.9   ...   ...   ...
...   ...   ...   ...
4.0   ...   ...   ...
...   ...   ...   ...
24.0  ...   ...   ...

Je veux ajuster les paramètres dans param, de sorte que les variables Y résultantes correspondent le mieux à mes valeurs de référence, par exemple :

t         Y(1)  Y(2)  Y(3)
0.5       1.1   N/A   N/A
1.0       1.9   N/A   N/A
4.0       2.3   2.7   2.1
5.0       N/A   2.6   2.2
24.0      0.9   1.5   2.0

Quelle routine Octave/Matlab (d'autres langages sont les bienvenus) peut effectuer un ajustement multiparamétrique (moindres carrés / spline) ? Comment est-il possible de combiner des ensembles de paramètres pour différentes valeurs initiales Y0 dans l'ajustement ? Je serais heureux si vous pouviez me donner des indications et des possibilités.

Je vous prie d'agréer, Madame, Monsieur, l'expression de mes salutations distinguées, Simon

Answer 1

Cela devrait être relativement simple avec scipy. scipy.optimize.leastsq() prend une fonction qui doit renvoyer un tableau de résidus pour un vecteur de paramètres donné. Elle minimisera la somme des carrés des résidus. Pour traiter plusieurs ensembles de données avec des valeurs initiales différentes, il suffit d'exécuter l'EDO une fois pour chaque ensemble de données, de calculer les résidus pour chaque paire ensemble de données/exécution, puis de concaténer les vecteurs de résidus ensemble. Voici un schéma approximatif :

import numpy
from scipy import integrate, optimize

# The initial guess.
p0 = numpy.array([7.2, 8.6, 9.5])

# The collected datasets.
# A list of (t, y0, y) tuples.
# The y's are all (len(y0), len(t))-shaped arrays. The output of
# integrate.odeint is also in this format.
datasets = [...]

def odes(y, t, params):
    dydt = [
        -params[0] * y[0],
        params[0]*y[0] - params[1]*y[1]*y[2],
        params[1]*y[1]*y[2],
    ]
    return np.array(dydt)

def residuals(params, datasets):
    res = []
    for t, y0, y in datasets:
        res.append(integrate.odeint(odes, y0, t, args=(params,)) - y)

    # Stack them horizontally and flatten the array into the expected vector.
    # You're on your own for handling missing data. Look into the numpy.ma
    # module.
    all_residuals = numpy.hstack(res).ravel()
    return all_residuals

opt_params, err = optimize.leastsq(residuals, p0, args=(datasets,))

Answer 2

Voulez-vous dire que chaque fonction y(t) doit être ajustée ? Dans ce cas, un ajustement par la méthode des moindres carrés ou par la méthode des splines pour chaque ensemble de Yi en fonction du temps fonctionnera parfaitement. Je ne peux pas dire lequel serait le meilleur sans voir vos données.

Vous devrez trouver une autre variable indépendante si vous voulez ajuster une courbe sur toutes les valeurs de Yi à un moment donné et observer l'évolution de cette courbe dans le temps.

MISE À JOUR : L'ajustement des moindres carrés est ce qu'il est - je n'ai pas de routine particulière à recommander. SciPy en a un Je n'en doute pas. Je suis désolé de ne pas avoir de meilleure recommandation. Je n'apprends que Python pour l'instant.

Je ne sais pas ce que vous entendez par "indicateur d'aptitude". L'ajustement des moindres carrés calcule les coefficients qui minimisent le carré moyen de l'erreur entre l'ajustement et les données en chaque point.

Pour combiner plusieurs ensembles de données en un seul ajustement, il suffit de les fusionner et de refaire le calcul.

Answer 3

J'ai développé une boîte à outils Matlab complète pour ajuster les paramètres et les valeurs initiales des ODE à de multiples ensembles de données expérimentales. Elle peut gérer différentes valeurs initiales en fonction de chaque expérience et est disponible à l'adresse suivante www.potterswheel.de .