Fitting multiple data with model and NDSolve with different initial conditions, and other shared parameters

8

5

I know that there are already questions about fitting multiple datasets and about NDSolve and about shared and non shared parameters, but I tried to apply them and some things are still not clear.

Here is my equation :

l = 10^(-5)
k = 1/l
chic = 0.5
T = 100

eq = {R'[t] == -a[t]*R[t] + b[t], 
  b'[t] == beta/2*(Tanh[(chi[t] - chic)*k] - 1), 
  a'[t] == -alpha/2*(Tanh[(chi[t] - chic)*k] - 1), 
  chi'[t] == -kappa*R[t]*(chi[t] - 2*chic), a[0] == a0, b[0] == b0, 
  R[0] == R0, chi[0] == 0}

I want to fit with regards to the variables : $alpha, beta, kappa, a0, b0$ as shared parameters and $R0$ as non shared parameter, meaning it would be different fr each one.

The joined data is given as an appendix just afterwards.

The non-joined data (meaning the 5 data-sets separately) looks like that :

enter image description here

So I tried to change $R0$ as a variable, and I got inspired by the answer of @JimB in Finding NonlinearModelFit of multiple data sets with the same parameters and in two dimensions :

model[alpha_?NumberQ, beta_?NumberQ, kappa_?NumberQ, a0_?NumberQ, 
  b0_?NumberQ] :=  (model[alpha, beta, kappa, a0, b0] = 
   Module[{R, chi, b, a, t, R0}, 
    First[R /. 
      NDSolve[{D[R[t, R0], t] == -a[t, R0]*R[t, R0] + 
          b[t, R0], 
        D[b[t, R0], t] == beta/2*(Tanh[(chi[t, R0] - chic)*k] - 1), 
        D[a[t, R0], t] == -alpha/2*(Tanh[(chi[t, R0] - chic)*k] - 1), 
        D[chi[t, R0], t] == -kappa*(chi[t, R0] - 2*chic), 
        a[0, R0] == a0, b[0, R0] == b0, R[0, R0] == R0, chi[0,R0] == 0}, {R, b, 
        a, chi}, {t, 0, T}, {R0, 0, 300}]]]);
nlm = NonlinearModelFit[data, 
   {model[alpha, beta, kappa, a0, b0][t, 
    R0], alpha >= 0, beta >= 0, kappa >= 0, a0 >= 0, b0 >= 0}, {{alpha, 0.1}, { beta, 0.1}, { kappa, 0.05}, {a0, 0.01}, {b0,
      3}}, {t, R0}];
nlm["BestFitParameters"]


The parameters are believed to be around :

alpha = 0.1
beta= 0.1
kappa = 0.05
a0 = 0.01
b0 = 3

But it didn't work... :

NonlinearModelFit::nrnum: The function value 1/2 ((-22.6124+R$3721[3.,22.])^2+(-119.51+R$3721[3.,119.])^2+(-24.738+R$3721[6.,22.])^2+(-60.1536+R$3721[6.,60.])^2+(-126.123+R$3721[6.,119.])^2+(-16.8895+R$3721[9.,17.])^2+(-25.4959+R$3721[9.,22.])^2+(-57.9807+R$3721[9.,60.])^2+(-110.446+R$3721[9.,119.])^2+(-17.3404+R$3721[12.,17.])^2+(-26.1946+R$3721[12.,22.])^2+(-60.9089+R$3721[12.,60.])^2+(-110.332+R$3721[12.,119.])^2+<<25>>+(-200.187+R$3721[27.,185.])^2+(-20.6519+R$3721[30.,17.])^2+(-34.5678+R$3721[30.,22.])^2+(-68.705+R$3721[30.,60.])^2+(-111.198+R$3721[30.,119.])^2+(-199.25+R$3721[30.,185.])^2+(-19.4591+R$3721[33.,17.])^2+(-35.9263+R$3721[33.,22.])^2+(-68.2107+R$3721[33.,60.])^2+(-109.903+R$3721[33.,119.])^2+(-198.411+R$3721[33.,185.])^2+(-20.6855+R$3721[36.,17.])^2+<<819>>) is not a real number at {alpha,beta,kappa,a0,b0} = {0.1,0.1,0.05,0.01,3.}.

I assume there is an issue with $R0$, but I don't get where.

How could I proceed ?

Also, I don't know how I could fix a priori the initial conditions for each fit in order to extract only the shared parameters.

DATA

MathematicaStackExchange doesn't give the possibility to enter to much characters. I can give only the joined data.

1. joined data with R0 as a variable

Here is the joined data.

data={{9., 17., 16.8895}, {12., 17., 17.3404}, {15., 17., 17.1633}, {18., 
  17., 19.3417}, {21., 17., 17.9899}, {24., 17., 19.9677}, {27., 17., 
  19.4362}, {30., 17., 20.6519}, {33., 17., 19.4591}, {36., 17., 
  20.6855}, {39., 17., 20.1952}, {42., 17., 21.9949}, {45., 17., 
  21.0234}, {48., 17., 22.7408}, {51., 17., 22.3908}, {54., 17., 
  25.0918}, {57., 17., 23.5989}, {60., 17., 26.0703}, {63., 17., 
  24.5605}, {66., 17., 27.2539}, {69., 17., 26.1619}, {72., 17., 
  28.4762}, {75., 17., 27.5854}, {78., 17., 29.8393}, {81., 17., 
  28.3553}, {84., 17., 30.3221}, {87., 17., 29.675}, {90., 17., 
  31.5653}, {93., 17., 30.5337}, {96., 17., 33.3734}, {99., 17., 
  31.6876}, {102., 17., 34.1503}, {105., 17., 33.3065}, {108., 17., 
  35.3291}, {111., 17., 33.9209}, {114., 17., 36.773}, {117., 17., 
  35.4094}, {120., 17., 41.5902}, {123., 17., 36.1305}, {126., 17., 
  37.971}, {129., 17., 36.402}, {132., 17., 39.1158}, {135., 17., 
  38.0177}, {138., 17., 40.8558}, {141., 17., 39.6065}, {144., 17., 
  40.9749}, {147., 17., 39.8896}, {150., 17., 41.8237}, {153., 17., 
  40.5802}, {156., 17., 42.3858}, {159., 17., 40.6619}, {162., 17., 
  44.4442}, {165., 17., 45.4162}, {168., 17., 46.1884}, {171., 17., 
  44.6008}, {174., 17., 47.1647}, {177., 17., 45.3808}, {180., 17., 
  46.5859}, {183., 17., 45.3035}, {186., 17., 47.6604}, {189., 17., 
  46.6771}, {192., 17., 45.9242}, {195., 17., 46.767}, {198., 17., 
  44.6899}, {201., 17., 46.6628}, {204., 17., 46.1571}, {207., 17., 
  46.5555}, {210., 17., 44.835}, {213., 17., 45.1423}, {216., 17., 
  45.1954}, {219., 17., 45.309}, {222., 17., 47.7791}, {225., 17., 
  46.7777}, {228., 17., 48.135}, {231., 17., 45.6493}, {234., 17., 
  45.8933}, {237., 17., 46.1803}, {240., 17., 46.7285}, {243., 17., 
  46.8063}, {246., 17., 47.1679}, {249., 17., 46.8787}, {252., 17., 
  47.2715}, {255., 17., 47.5362}, {258., 17., 48.9234}, {261., 17., 
  47.5456}, {264., 17., 53.5554}, {267., 17., 52.5704}, {270., 17., 
  49.6049}, {273., 17., 49.1189}, {276., 17., 48.9498}, {279., 17., 
  49.6024}, {282., 17., 49.7491}, {285., 17., 53.1681}, {288., 17., 
  51.7124}, {291., 17., 50.8069}, {294., 17., 50.0237}, {297., 17., 
  50.5922}, {300., 17., 50.6518}, {303., 17., 50.8827}, {306., 17., 
  51.2245}, {309., 17., 51.0911}, {312., 17., 52.3379}, {315., 17., 
  52.5112}, {318., 17., 53.9182}, {321., 17., 53.7082}, {324., 17., 
  54.9239}, {327., 17., 53.7369}, {330., 17., 51.7204}, {333., 17., 
  55.993}, {336., 17., 56.8489}, {339., 17., 53.3037}, {342., 17., 
  52.0201}, {345., 17., 52.6267}, {348., 17., 52.5615}, {351., 17., 
  55.4133}, {354., 17., 55.5549}, {357., 17., 52.2672}, {360., 17., 
  54.2202}, {363., 17., 50.3245}, {366., 17., 54.0435}, {369., 17., 
  51.0724}, {372., 17., 51.2091}, {375., 17., 51.6602}, {378., 17., 
  51.3684}, {381., 17., 51.5346}, {384., 17., 51.9204}, {387., 17., 
  52.3952}, {390., 17., 52.9114}, {393., 17., 54.3833}, {396., 17., 
  55.1898}, {399., 17., 51.3853}, {402., 17., 55.048}, {405., 17., 
  50.8574}, {408., 17., 51.9619}, {411., 17., 52.5775}, {414., 17., 
  52.5676}, {417., 17., 51.0891}, {420., 17., 54.3895}, {423., 17., 
  54.7591}, {426., 17., 53.9934}, {429., 17., 53.8877}, {435., 17., 
  55.4067}, {441., 17., 56.0656}, {447., 17., 57.4607}, {453., 17., 
  51.6628}, {456., 17., 54.3568}, {459., 17., 57.6827}, {465., 17., 
  54.8474}, {468., 17., 51.0797}, {471., 17., 53.1862}, {474., 17., 
  53.3921}, {477., 17., 54.468}, {480., 17., 54.1083}, {483., 17., 
  50.7948}, {486., 17., 53.3431}, {489., 17., 48.8646}, {492., 17., 
  53.3906}, {495., 17., 51.6016}, {498., 17., 54.1742}, {501., 17., 
  54.6549}, {504., 17., 50.0598}, {507., 17., 53.849}, {510., 17., 
  52.6431}, {513., 17., 54.3103}, {516., 17., 50.5004}, {519., 17., 
  50.8213}, {522., 17., 50.8512}, {525., 17., 52.4319}, {528., 17., 
  55.2716}, {3., 22., 22.6124}, {6., 22., 24.738}, {9., 22., 
  25.4959}, {12., 22., 26.1946}, {15., 22., 27.6091}, {18., 22., 
  29.1024}, {21., 22., 30.6462}, {24., 22., 32.9126}, {27., 22., 
  34.1471}, {30., 22., 34.5678}, {33., 22., 35.9263}, {36., 22., 
  37.4284}, {39., 22., 38.5027}, {42., 22., 39.5611}, {45., 22., 
  40.743}, {48., 22., 41.9482}, {51., 22., 42.7558}, {54., 22., 
  43.5064}, {57., 22., 44.43}, {60., 22., 45.7449}, {63., 22., 
  47.0524}, {66., 22., 48.0848}, {69., 22., 48.8836}, {72., 22., 
  49.6807}, {75., 22., 50.6801}, {78., 22., 51.6959}, {81., 22., 
  52.6475}, {84., 22., 53.5902}, {87., 22., 54.4008}, {90., 22., 
  54.774}, {93., 22., 55.6085}, {96., 22., 56.3299}, {99., 22., 
  56.4428}, {102., 22., 56.7936}, {105., 22., 57.4926}, {108., 22., 
  58.2406}, {111., 22., 59.1169}, {114., 22., 59.5766}, {117., 22., 
  59.7909}, {120., 22., 61.6917}, {123., 22., 62.4342}, {126., 22., 
  61.5979}, {129., 22., 61.8203}, {132., 22., 62.5629}, {135., 22., 
  63.4556}, {138., 22., 63.688}, {141., 22., 63.9159}, {144., 22., 
  63.9802}, {147., 22., 64.1833}, {150., 22., 64.3304}, {153., 22., 
  64.3847}, {156., 22., 64.6173}, {159., 22., 64.9009}, {162., 22., 
  65.1622}, {165., 22., 65.4684}, {168., 22., 65.5182}, {171., 22., 
  66.1171}, {174., 22., 66.4103}, {177., 22., 66.2592}, {180., 22., 
  66.185}, {183., 22., 65.8147}, {186., 22., 65.733}, {189., 22., 
  65.6618}, {192., 22., 64.7882}, {195., 22., 64.8274}, {198., 22., 
  64.9444}, {201., 22., 63.1305}, {204., 22., 62.3995}, {207., 22., 
  63.0431}, {210., 22., 62.2181}, {213., 22., 62.5286}, {216., 22., 
  62.1711}, {219., 22., 60.8353}, {222., 22., 60.7586}, {225., 22., 
  60.7004}, {228., 22., 59.5638}, {231., 22., 59.1517}, {234., 22., 
  58.9346}, {237., 22., 59.0493}, {240., 22., 59.5229}, {243., 22., 
  58.0876}, {246., 22., 56.247}, {249., 22., 56.173}, {252., 22., 
  56.1419}, {255., 22., 55.2417}, {258., 22., 56.2456}, {261., 22., 
  57.9169}, {264., 22., 60.728}, {267., 22., 63.6912}, {270., 22., 
  61.9647}, {273., 22., 57.0852}, {276., 22., 54.2803}, {279., 22., 
  55.3487}, {282., 22., 58.0208}, {285., 22., 60.8749}, {288., 22., 
  61.029}, {291., 22., 59.3053}, {294., 22., 56.7078}, {297., 22., 
  53.8873}, {300., 22., 55.2545}, {303., 22., 56.5482}, {306., 22., 
  56.0664}, {309., 22., 55.2537}, {312., 22., 55.3196}, {315., 22., 
  55.8909}, {318., 22., 55.6318}, {321., 22., 56.213}, {324., 22., 
  55.4207}, {327., 22., 54.2877}, {330., 22., 55.1178}, {333., 22., 
  51.193}, {336., 22., 48.5713}, {339., 22., 49.5028}, {342., 22., 
  49.4166}, {345., 22., 50.0304}, {348., 22., 50.9326}, {351., 22., 
  52.014}, {354., 22., 50.2956}, {357., 22., 49.8529}, {360., 22., 
  50.8205}, {363., 22., 51.376}, {366., 22., 50.6679}, {369., 22., 
  51.6815}, {372., 22., 53.5813}, {375., 22., 53.7359}, {378., 22., 
  54.6252}, {381., 22., 55.2786}, {384., 22., 53.4308}, {387., 22., 
  54.5401}, {390., 22., 57.9795}, {393., 22., 55.2026}, {396., 22., 
  55.386}, {399., 22., 59.8766}, {402., 22., 58.1028}, {405., 22., 
  57.129}, {408., 22., 56.9853}, {411., 22., 57.2221}, {414., 22., 
  56.9648}, {417., 22., 55.586}, {420., 22., 56.7903}, {423., 22., 
  56.2825}, {426., 22., 53.8012}, {429., 22., 52.6652}, {432., 22., 
  54.2455}, {435., 22., 56.3002}, {438., 22., 56.2343}, {441., 22., 
  56.7575}, {444., 22., 56.7977}, {447., 22., 56.3049}, {450., 22., 
  54.6538}, {453., 22., 52.5136}, {456., 22., 52.3433}, {459., 22., 
  52.828}, {462., 22., 54.0433}, {465., 22., 51.5131}, {468., 22., 
  50.4781}, {471., 22., 52.6831}, {474., 22., 52.4475}, {477., 22., 
  52.6825}, {480., 22., 52.5579}, {483., 22., 52.8213}, {486., 22., 
  53.6997}, {489., 22., 53.3714}, {492., 22., 52.3218}, {495., 22., 
  52.3176}, {498., 22., 53.8036}, {501., 22., 53.7502}, {504., 22., 
  55.6969}, {507., 22., 56.1864}, {510., 22., 52.9824}, {513., 22., 
  55.2477}, {516., 22., 54.727}, {519., 22., 54.0447}, {522., 22., 
  56.1034}, {525., 22., 53.0694}, {528., 22., 51.3001}, {6., 60., 
  60.1536}, {9., 60., 57.9807}, {12., 60., 60.9089}, {15., 60., 
  59.4291}, {18., 60., 61.3227}, {21., 60., 61.8788}, {24., 60., 
  67.2192}, {27., 60., 66.2767}, {30., 60., 68.705}, {33., 60., 
  68.2107}, {36., 60., 70.8731}, {39., 60., 68.7269}, {42., 60., 
  73.2306}, {45., 60., 72.3068}, {48., 60., 74.8006}, {51., 60., 
  72.1975}, {54., 60., 76.577}, {57., 60., 75.5894}, {60., 60., 
  76.342}, {63., 60., 75.5134}, {66., 60., 77.47}, {69., 60., 
  76.6854}, {72., 60., 78.7422}, {75., 60., 78.6074}, {78., 60., 
  81.0158}, {81., 60., 82.8521}, {84., 60., 85.1395}, {87., 60., 
  85.211}, {90., 60., 84.5157}, {93., 60., 83.622}, {96., 60., 
  88.1703}, {99., 60., 85.6195}, {102., 60., 86.8345}, {105., 60., 
  86.5568}, {108., 60., 87.5942}, {111., 60., 88.3053}, {114., 60., 
  88.3475}, {117., 60., 89.3993}, {120., 60., 91.7091}, {123., 60., 
  89.7268}, {126., 60., 90.6704}, {129., 60., 89.7999}, {132., 60., 
  90.369}, {135., 60., 88.7787}, {138., 60., 90.3022}, {141., 60., 
  89.8267}, {144., 60., 91.2241}, {147., 60., 91.2859}, {150., 60., 
  92.992}, {153., 60., 91.0079}, {156., 60., 93.0784}, {159., 60., 
  90.8868}, {162., 60., 92.7426}, {165., 60., 92.757}, {168., 60., 
  94.4202}, {171., 60., 92.2914}, {174., 60., 90.3876}, {177., 60., 
  89.3376}, {180., 60., 89.814}, {183., 60., 88.9134}, {186., 60., 
  89.7058}, {189., 60., 91.642}, {192., 60., 90.3205}, {195., 60., 
  87.8566}, {198., 60., 87.6065}, {201., 60., 87.0403}, {204., 60., 
  87.3344}, {207., 60., 87.2313}, {210., 60., 87.3705}, {213., 60., 
  86.9135}, {216., 60., 87.2684}, {219., 60., 87.2989}, {222., 60., 
  85.4766}, {225., 60., 85.3534}, {228., 60., 86.535}, {231., 60., 
  86.1929}, {234., 60., 86.089}, {237., 60., 85.9466}, {240., 60., 
  85.1389}, {243., 60., 85.0242}, {246., 60., 84.4313}, {249., 60., 
  83.7604}, {252., 60., 81.9419}, {255., 60., 83.773}, {258., 60., 
  82.7046}, {261., 60., 84.7331}, {264., 60., 86.0393}, {267., 60., 
  84.7472}, {270., 60., 79.1677}, {273., 60., 80.9426}, {276., 60., 
  79.9624}, {279., 60., 75.5272}, {282., 60., 79.3103}, {285., 60., 
  80.8015}, {288., 60., 81.3927}, {291., 60., 80.1678}, {294., 60., 
  80.268}, {297., 60., 79.9067}, {300., 60., 76.9766}, {303., 60., 
  81.8132}, {306., 60., 73.6449}, {309., 60., 76.4059}, {312., 60., 
  76.4056}, {315., 60., 81.7311}, {318., 60., 80.8468}, {321., 60., 
  80.958}, {324., 60., 86.9248}, {327., 60., 78.3434}, {330., 60., 
  74.8752}, {333., 60., 78.0912}, {336., 60., 81.5165}, {339., 60., 
  72.7919}, {342., 60., 74.2966}, {345., 60., 79.2233}, {348., 60., 
  81.9791}, {351., 60., 74.3276}, {354., 60., 85.1221}, {357., 60., 
  78.6944}, {360., 60., 75.8183}, {363., 60., 75.6696}, {366., 60., 
  75.9147}, {369., 60., 76.3326}, {372., 60., 80.0048}, {375., 60., 
  79.8311}, {378., 60., 79.0427}, {381., 60., 81.8084}, {384., 60., 
  73.5742}, {387., 60., 84.2291}, {390., 60., 84.9122}, {393., 60., 
  82.6657}, {396., 60., 78.2888}, {399., 60., 90.0235}, {402., 60., 
  83.3667}, {405., 60., 81.7737}, {408., 60., 81.19}, {411., 60., 
  82.3131}, {414., 60., 79.8072}, {417., 60., 74.4822}, {420., 60., 
  75.6291}, {423., 60., 82.2655}, {426., 60., 73.704}, {429., 60., 
  81.4184}, {432., 60., 72.1127}, {435., 60., 74.7053}, {438., 60., 
  79.4664}, {441., 60., 86.4491}, {444., 60., 79.5096}, {447., 60., 
  77.1761}, {450., 60., 83.082}, {453., 60., 80.3418}, {456., 60., 
  85.3873}, {459., 60., 85.7409}, {462., 60., 73.3735}, {465., 60., 
  72.2276}, {468., 60., 82.7752}, {471., 60., 71.6917}, {474., 60., 
  78.5233}, {477., 60., 82.4042}, {480., 60., 83.8073}, {483., 60., 
  91.5845}, {486., 60., 82.8906}, {489., 60., 87.3935}, {492., 60., 
  89.9856}, {495., 60., 74.1819}, {498., 60., 77.5752}, {501., 60., 
  82.6796}, {504., 60., 79.2659}, {507., 60., 81.5865}, {510., 60., 
  82.709}, {513., 60., 88.4083}, {516., 60., 81.7317}, {519., 60., 
  76.2638}, {522., 60., 86.2863}, {525., 60., 93.2163}, {528., 60., 
  82.6943}, {3., 119., 119.51}, {6., 119., 126.123}, {9., 119., 
  110.446}, {12., 119., 110.332}, {15., 119., 110.478}, {18., 119., 
  111.335}, {21., 119., 109.536}, {24., 119., 111.901}, {27., 119., 
  110.46}, {30., 119., 111.198}, {33., 119., 109.903}, {36., 119., 
  110.72}, {39., 119., 109.635}, {42., 119., 110.643}, {45., 119., 
  109.528}, {48., 119., 110.348}, {51., 119., 110.117}, {54., 119., 
  109.117}, {57., 119., 108.536}, {60., 119., 108.615}, {63., 119., 
  109.495}, {66., 119., 111.304}, {69., 119., 111.139}, {72., 119., 
  114.285}, {75., 119., 113.627}, {78., 119., 114.77}, {81., 119., 
  114.544}, {84., 119., 115.304}, {87., 119., 114.895}, {90., 119., 
  115.859}, {93., 119., 114.357}, {96., 119., 115.038}, {99., 119., 
  114.305}, {102., 119., 115.09}, {105., 119., 114.815}, {108., 119., 
  113.203}, {111., 119., 113.46}, {114., 119., 114.573}, {117., 119., 
  113.339}, {120., 119., 114.354}, {123., 119., 112.285}, {126., 119.,
   112.695}, {129., 119., 112.032}, {132., 119., 112.253}, {135., 
  119., 108.945}, {138., 119., 109.271}, {141., 119., 108.654}, {144.,
   119., 104.336}, {147., 119., 103.609}, {150., 119., 
  105.778}, {153., 119., 105.077}, {156., 119., 104.868}, {159., 119.,
   103.945}, {162., 119., 104.039}, {165., 119., 101.727}, {168., 
  119., 97.6562}, {171., 119., 99.6703}, {174., 119., 96.6503}, {177.,
   119., 98.3032}, {180., 119., 98.8859}, {183., 119., 
  97.9825}, {186., 119., 94.8383}, {189., 119., 93.4101}, {192., 119.,
   88.9132}, {195., 119., 91.7409}, {198., 119., 93.2425}, {201., 
  119., 86.1268}, {204., 119., 84.9263}, {207., 119., 86.3445}, {210.,
   119., 84.4667}, {213., 119., 85.9353}, {216., 119., 
  85.7998}, {219., 119., 85.2672}, {222., 119., 86.3356}, {225., 119.,
   86.7423}, {228., 119., 86.1353}, {231., 119., 84.8631}, {234., 
  119., 84.7305}, {237., 119., 83.385}, {240., 119., 87.5174}, {243., 
  119., 83.3014}, {246., 119., 86.9219}, {249., 119., 78.3219}, {252.,
   119., 78.9197}, {255., 119., 74.785}, {258., 119., 67.8261}, {261.,
   119., 75.8036}, {264., 119., 86.2339}, {267., 119., 
  87.3689}, {270., 119., 88.1322}, {273., 119., 86.1332}, {276., 119.,
   89.9111}, {279., 119., 90.5619}, {282., 119., 88.4012}, {285., 
  119., 85.5809}, {288., 119., 76.692}, {291., 119., 80.0753}, {294., 
  119., 90.1118}, {297., 119., 91.8565}, {300., 119., 85.0882}, {303.,
   119., 89.1269}, {306., 119., 96.8869}, {309., 119., 
  75.4618}, {312., 119., 96.3013}, {315., 119., 89.4435}, {318., 119.,
   103.21}, {321., 119., 94.6233}, {324., 119., 102.48}, {327., 119., 
  96.7664}, {330., 119., 84.2408}, {333., 119., 97.3822}, {336., 119.,
   74.2619}, {339., 119., 87.2886}, {342., 119., 118.024}, {345., 
  119., 113.648}, {348., 119., 112.4}, {351., 119., 107.295}, {354., 
  119., 111.618}, {357., 119., 112.181}, {360., 119., 112.119}, {363.,
   119., 90.6252}, {366., 119., 106.837}, {369., 119., 
  99.7227}, {372., 119., 97.5255}, {375., 119., 108.211}, {378., 119.,
   117.211}, {381., 119., 97.9301}, {384., 119., 104.567}, {387., 
  119., 117.343}, {390., 119., 121.622}, {393., 119., 106.117}, {396.,
   119., 116.022}, {399., 119., 118.856}, {402., 119., 
  106.854}, {405., 119., 112.418}, {408., 119., 112.79}, {411., 119., 
  112.225}, {414., 119., 116.686}, {417., 119., 111.297}, {420., 119.,
   115.404}, {423., 119., 117.563}, {426., 119., 116.243}, {429., 
  119., 119.805}, {432., 119., 112.863}, {435., 119., 103.505}, {438.,
   119., 116.846}, {441., 119., 115.508}, {444., 119., 
  115.579}, {447., 119., 101.756}, {450., 119., 102.848}, {453., 119.,
   112.506}, {456., 119., 113.93}, {459., 119., 116.386}, {462., 119.,
   108.138}, {465., 119., 108.635}, {468., 119., 110.514}, {471., 
  119., 108.217}, {474., 119., 110.008}, {477., 119., 95.7788}, {480.,
   119., 92.8073}, {483., 119., 104.382}, {486., 119., 98.77}, {489., 
  119., 112.527}, {492., 119., 94.6092}, {495., 119., 89.2861}, {498.,
   119., 92.0002}, {501., 119., 98.7618}, {504., 119., 
  105.274}, {507., 119., 96.7057}, {510., 119., 93.5207}, {513., 119.,
   90.5992}, {516., 119., 87.1486}, {519., 119., 103.466}, {522., 
  119., 100.133}, {525., 119., 120.605}, {528., 119., 125.717}, {12., 
  185., 185.791}, {15., 185., 199.035}, {18., 185., 197.796}, {21., 
  185., 185.256}, {24., 185., 199.576}, {27., 185., 200.187}, {30., 
  185., 199.25}, {33., 185., 198.411}, {36., 185., 198.288}, {39., 
  185., 194.506}, {42., 185., 189.658}, {45., 185., 191.203}, {48., 
  185., 185.757}, {51., 185., 183.642}, {54., 185., 183.513}, {57., 
  185., 186.524}, {60., 185., 182.793}, {63., 185., 182.218}, {66., 
  185., 182.045}, {69., 185., 176.614}, {72., 185., 182.432}, {75., 
  185., 181.409}, {78., 185., 182.438}, {81., 185., 179.939}, {84., 
  185., 182.435}, {87., 185., 181.521}, {90., 185., 176.654}, {93., 
  185., 175.39}, {96., 185., 179.446}, {99., 185., 173.541}, {102., 
  185., 176.645}, {105., 185., 176.715}, {108., 185., 173.915}, {111.,
   185., 173.14}, {114., 185., 173.045}, {117., 185., 160.089}, {120.,
   185., 165.306}, {123., 185., 165.906}, {126., 185., 
  165.712}, {129., 185., 159.285}, {132., 185., 163.219}, {135., 185.,
   156.287}, {138., 185., 150.445}, {141., 185., 153.388}, {144., 
  185., 138.083}, {147., 185., 137.152}, {150., 185., 133.003}, {153.,
   185., 130.634}, {156., 185., 131.832}, {159., 185., 
  136.142}, {162., 185., 133.906}, {165., 185., 130.929}, {168., 185.,
   136.717}, {171., 185., 129.749}, {174., 185., 148.377}, {177., 
  185., 133.068}, {180., 185., 149.921}, {183., 185., 134.802}, {186.,
   185., 150.543}, {189., 185., 138.678}, {192., 185., 147.06}, {195.,
   185., 143.604}, {198., 185., 143.368}, {201., 185., 
  140.587}, {204., 185., 138.171}, {207., 185., 140.699}, {210., 185.,
   137.346}, {213., 185., 126.241}, {216., 185., 131.743}, {219., 
  185., 134.835}, {222., 185., 134.086}, {225., 185., 137.185}, {228.,
   185., 135.892}, {231., 185., 141.62}, {234., 185., 135.963}, {237.,
   185., 133.382}, {240., 185., 134.258}, {243., 185., 
  141.568}, {246., 185., 137.642}, {249., 185., 131.681}, {252., 185.,
   132.635}, {255., 185., 134.506}, {258., 185., 136.089}, {261., 
  185., 138.973}, {264., 185., 141.048}, {267., 185., 133.785}, {270.,
   185., 133.245}, {273., 185., 116.408}, {276., 185., 123.9}, {279., 
  185., 120.251}, {282., 185., 116.984}, {285., 185., 135.753}, {288.,
   185., 123.026}, {291., 185., 112.116}, {294., 185., 
  134.164}, {297., 185., 134.548}, {300., 185., 129.032}, {303., 185.,
   116.97}, {306., 185., 113.993}, {309., 185., 99.4695}, {312., 185.,
   97.4854}, {315., 185., 100.422}, {318., 185., 117.461}, {321., 
  185., 99.4758}, {324., 185., 106.366}, {327., 185., 108.271}, {330.,
   185., 104.738}, {333., 185., 117.487}, {336., 185., 
  101.704}, {339., 185., 101.32}, {342., 185., 112.97}, {345., 185., 
  96.6092}, {348., 185., 99.2531}, {351., 185., 120.19}, {354., 185., 
  124.284}, {357., 185., 130.082}, {360., 185., 121.699}, {363., 185.,
   108.539}, {366., 185., 103.98}, {369., 185., 100.293}, {372., 185.,
   94.7848}, {375., 185., 103.281}, {378., 185., 114.4}, {381., 185., 
  94.8752}, {384., 185., 101.51}, {387., 185., 104.285}, {390., 185., 
  107.424}, {393., 185., 112.506}, {396., 185., 104.061}, {399., 185.,
   113.713}, {402., 185., 136.378}, {405., 185., 134.92}, {408., 185.,
   139.111}, {411., 185., 143.397}, {414., 185., 139.998}, {417., 
  185., 137.19}, {420., 185., 143.812}, {423., 185., 133.346}, {426., 
  185., 141.8}, {429., 185., 136.171}, {432., 185., 137.842}, {435., 
  185., 147.509}, {438., 185., 140.488}, {441., 185., 142.855}, {444.,
   185., 151.992}, {447., 185., 145.348}, {450., 185., 
  138.757}, {453., 185., 135.964}, {456., 185., 140.381}, {459., 185.,
   143.697}, {462., 185., 136.854}, {465., 185., 129.477}, {468., 
  185., 138.181}, {471., 185., 142.726}, {474., 185., 143.633}, {477.,
   185., 133.913}, {480., 185., 157.635}, {483., 185., 
  147.941}, {486., 185., 142.015}, {489., 185., 130.545}, {492., 185.,
   141.941}, {495., 185., 142.863}, {498., 185., 135.462}, {501., 
  185., 139.637}, {504., 185., 128.002}, {507., 185., 140.211}, {510.,
   185., 140.209}, {513., 185., 132.36}, {516., 185., 141.088}, {519.,
   185., 142.756}, {522., 185., 152.256}, {525., 185., 
  164.725}, {528., 185., 153.737}}

J.A

Posted 2019-05-24T10:05:34.487

Reputation: 1 193

1Observe the output of the following code piece and think about the reason: expr = x; Module[{x}, {x, expr}] – xzczd – 2019-05-26T11:36:20.830

@xzczd Thanks for helping. Do you mean that x is going to be evaluated and expr is not going to be evaluated and that's why the output in your example is {x$24532, x} ? – J.A – 2019-05-26T12:24:52.853

I have to admit I do not understand the link between this and my issue.... Could you be a bit more explicit plz ? – J.A – 2019-05-26T12:32:19.630

1No, this isn't that related to evaluation order. (Though it's not impossible to fix by tackling evaluation order. ) The main issue here is, 1. Module localizes variable by renaming it as …$…. (In the case above, x becomes x$24532); 2. When localizing, Module only sees explicit variable i.e. the x inside expr is not noticed by Module so it's not localized and not renamed. And the same thing happens on your R. – xzczd – 2019-05-26T12:39:20.117

1Only those explicit Rs are localized. You have Rs inside eq. – xzczd – 2019-05-26T12:48:20.363

So I replaced eq by an explicit formulation, and I entered chi,a,b into the local variables specification, but it's still not working. – J.A – 2019-05-26T12:54:00.503

1Please show the specific code. Also, the {R0, 0, 300} in NDSolve only makes the code unnecessarily complicated, and you've missed the initial condition for chi in your code, and D[chi[t, R0], t][t, R0] is obviously wrong. As to the update, alpha >= 0, beta >= 0, kappa >= 0, a0 >= 0, b0 >= 0 doesn't make sense, please don't guess the syntax, check the document of NDSolve carefully. – xzczd – 2019-05-26T13:17:07.803

Let us continue this discussion in chat.

– J.A – 2019-05-26T13:49:23.600

If I ListPlot your data I do not get your figure. In fact, data is a list of $3$-tuples, so I'm not sure how you expect to reconstruct the five graphs in your figure. Can you clarify how to obtain the figure from data? What is the exact code you use to generate it? – AccidentalFourierTransform – 2019-05-26T20:48:20.960

Are we supposed to fit the joined data or the non-joined data? If the latter, can you please post it? – AccidentalFourierTransform – 2019-05-26T21:26:00.970

@AccidentalFourierTransform I'm sorry but there would be too much characters, and the stackexchange doesn't allow it... The method I used used the joined data, with R0 as a variable. – J.A – 2019-05-26T21:48:50.943

1@J.A You can (and arguably should) use pastebin.com. – AccidentalFourierTransform – 2019-05-26T22:18:46.703

1

@AccidentalFourierTransform Another possible choice is to use the selected cell button of SE uploader.

– xzczd – 2019-05-27T05:25:17.287

Answers

7

The subject of parameter fitting comes up frequently on MSE. Parameter fitting is a difficult subject and will depend on your data quality, your model, and your intial guesses. I have been dabbling with StringTemplates as a potential way to encapsulate some of the basic parameter fitting work flow.

Approach

  • Use ParametricNDSolveValue to create the model.
  • Use StringTemplates to handle lists of parameters and variables.
  • Generate a Manipulate slider model to debug model and understand the effects of parameter changes.
  • Transfer initial guesses from manipulate to perform a fit.

Implementation

I commented the code so I hope it self explanatory. First assign the constants and prep the data.

(* Evaluate data first *)
(* Constants *)
l = 10^(-5);
k = 1/l;
chic = 0.5;
T = 550;
(* Get unique R0s *)
R0s = Union@data[[All, 2]];
(* Subset Matching R0 and Delete 2nd Column *)
rdat = (Cases[data, {_, #, _}][[All, {1, 3}]] & /@ R0s);

Now, set up the equations and the Manipulate slider to view how the model behaves and try to improve initial parameters estimates.

(* Generate System of Differential Equations *)
e1 = R'[t] == -a[t]*R[t] + b[t];
e3 = b'[t] == beta/2*(Tanh[(chi[t] - chic)*k] - 1);
e2 = a'[t] == -alpha/2*(Tanh[(chi[t] - chic)*k] - 1);
e4 = chi'[t] == -kappa*R[t]*(chi[t] - 2*chic);
ics = {a[0] == a0, b[0] == b0, R[0] == R0, chi[0] == 0};
eqns = {e1, e2, e3, e4}~Join~ics;
(*Variables*)
vbles = {R, a, b, chi};
(*Parameters with target and desired ranges*)
mat = {
   {alpha, 0.1, 0.00025, 0.5},
   {beta, 0.1, 0.00025, 0.5},
   {kappa, 0.05, 0.0125, 0.1},
   {a0, 0.01, 0.00005, 0.1},
   {b0, 3, 1, 6},
   {R0, 17, 17, 185}
   };
(* reduce the matrix because R0 does not participate in parameter \
fits *)
rmat = mat[[1 ;; -2]];
(* Build Manipulate sliders *)
sfun =  StringRiffle[(StringTemplate[
         "{{`1`,`2`},`3`,`4`,Appearance\[Rule]\"Labeled\"}"] @@ #) & \
/@ #, ","] &;
sliders = sfun[rmat];
(* Extract Parameters from mat *)
parms = mat[[All, 1]];
rparms = rmat[[All, 1]];
(* Create String Representations of parms *)
sparms = StringRiffle[ToString[#] & /@ parms, ","];
rsparms = StringRiffle[ToString[#] & /@ rparms, ","];
(* Create patterns and string reps of parameters *)
pats = Pattern @@@ (#*_ & /@ parms);
spats = StringRiffle[ToString[#] & /@ pats, ","];
(* List Plot of the data *)
lp = Graphics[{Hue[#2/185], PointSize[0.01], Point[{#1, #3}]} & @@@ 
    data, Axes -> True];
(* ParametricNDSolveValue *)
pfun = ParametricNDSolveValue[eqns, vbles, {t, 0, T}, parms];
(*Create an appropriate model function to fit*)
modelstring = "(#[[1]])&";
(* Create some PlotLegends *)
pl = ",PlotLegends\[Rule]{" <> 
   StringRiffle["\"R0=" <> ToString[#] <> "\"" & /@ R0s, ","] <> "}";
(* Build the model expression *)
ToExpression[
 StringTemplate[
   "model[`pats`][t_]:=`ms`@Through[pfun[`params`][t],List]\
/;And@@NumericQ/@{`params`};"][<|"pats" -> spats, "params" -> sparms, 
   "ms" -> modelstring|>]]
(* Create slider model *)
globalstring = 
  StringTemplate["global={`params`};"][<|"params" -> rsparms|>];
mantemp = 
  "Manipulate[`g`\[IndentingNewLine]Show[lp,Plot[Evaluate@({model[\
alpha,beta,kappa,a0,b0,#][t]}&/@R0s),{t,0,T},PlotRange\[Rule]{0,200}`\
pl`],ImageSize->Large],`sliders`]";
ToExpression@
 StringTemplate[mantemp][<|"sliders" -> sliders, "params" -> rsparms, 
   "pl" -> pl, "g" -> globalstring|>]
(*Display global variable*)
Dynamic@global

Manipulate Image

Now set up to fit the funtions for each R0 value.

(* Grab The initial parameter guesses *)
initguess = MapThread[List, {rparms, First@Dynamic@global}];
(* Create a fit function to operate on different R0s *)
fitfn = FindFit[rdat[[#]], 
    model[alpha, beta, kappa, a0, b0, R0s[[#]]][t], initguess, t, 
    Method -> "Gradient"] &;
(* Perform Fits on R0s *)
fits = fitfn[#][[All, 2]] & /@ Range@Length@R0s;
(* Display Results *)
fits // MatrixForm
Mean@fits

Parameter Estimates

The data is noisy leading to some dodgy results for the high R0. You can experiment with different fitting options, but you may need to improve your model and/or your data acquisition.

Update to Fit per $R_0$ Data set

As requested, here is a way to fit per data set. I also allowed $R_0$ to be fit using the column value as an initial guess. In this case, each fitted row is plotted. A word of caution, some fitting methods will run forever, so you may need to experiment.

(* Grab The initial parameter guesses from dynamic variable of slider \
*)
initguess = 
  MapThread[List, {parms, (First@Dynamic@global)~Join~{R0s[[#]]}}] &;
(* Create a fit function to operate on different R0s *)
fitfn = FindFit[rdat[[#]], model[alpha, beta, kappa, a0, b0, R0][t], 
    initguess[#], t, Method -> "Gradient", WorkingPrecision -> 10] &;
(* Perform Fits on R0s *)
(*fits = fitfn[#][[All,2]]&/@Range@Length@R0s;*)
fits = fitfn[#][[All, 2]] & /@ {1, 2, 3, 4, 5};
(* Display Results *)
fits // MatrixForm
mfit = Mean@fits
mat2 = rmat;
mat2[[All, 2]] = mfit[[1 ;; -2]];
Show[{lp, 
  Plot[Evaluate@((model @@ #)[t] & /@ fits), {t, 0, T}, 
   PlotRange -> {0, 200}, 
   PlotLegends -> {"R0=17.", "R0=22.", "R0=60.", "R0=119.", 
     "R0=185."}]}, ImageSize -> Large]

Individual Fits

Tim Laska

Posted 2019-05-24T10:05:34.487

Reputation: 12 026

Thanks this is a great job, and a very nice way to proceed. The only question I have is with regard to the fitting process. Is taking the mean the only way to proceed ? Meaning is there another way to gather the informations from the different data-sets ? – J.A – 2019-05-27T15:43:11.340

I added a plot of the fitted parameters for each row. – Tim Laska – 2019-05-27T21:06:13.807