1

I am trying to simultaneously fit two data sets that are separate functions of $x$ and $t$ with NonlinearModelFit.

Here is an example of the data:

datax = {{1, 2}, {2, 3}, {3, 4}, {4, 5}};
datat = {{1, 4}, {2, 7}, {3, 10}, {4, 13}};

Here is an example of models:

modelx[a_,b_,x_]:=a*x+b;
modelt[a_,b_,t_]:=a*t+b;

The parameters a and b are shared across the two models.

Here is an example of the code I'd hoped to use:

ResourceFunction["MultiNonlinearModelFit"][Join[{datax, datat}],
 {modelx[a, b, x], modelt[a, b, t]}, {a, b}, {x, t}];

The resulting error message is:

NonlinearModelFit: Number of coordinates (2) is not equal to the number of variables (3).

Edit

Consider the two data sets below. The first is a decaying exponential, the second is a polynomial. The first is a function of time, the second is a function of position.

nPoints = 256;
noiseLevel = 10^-1.5;

dataDecay = Table[{t, Sin[2 [Pi] 1.0 t] Exp[-0.5 t] + RandomReal[NormalDistribution[0, noiseLevel]]}, {t, 0, 10, 10/nPoints}] // Chop;

dataPolys = Table[{x, 1 + (1.0 x^2 + 0.5 x)/100 + RandomReal[NormalDistribution[0, noiseLevel]]}, {x, 0, 10, 10/nPoints}] // Chop;

I would like to create data sets that indicate whether the data belong to time or position:

dataDec = ArrayFlatten[{{Transpose[{ConstantArray[1, Length[dataDecay]]}], dataDecay}}];
dataPol = ArrayFlatten[{{Transpose[{ConstantArray[2, Length[dataPolys]]}], dataPolys}}];

I then model the data with the functions that created them:

modelDec[a_, b_, t_] := Sin[2 \[Pi] a t] Exp[-b t];
modelPol[a_, b_, x_] := 1 + (a x^2 + b x)/100;

So, in theory, when fitting the data, a should come out near 1 and b should come out near 0.5.

fit = NonlinearModelFit[Join[dataDec, dataPol], {modelDec[a, b, t] + modelPol[a, b, x]}, {{a, 1.0}, {b, 0.5}}, {t, x}, MaxIterations -> Infinity];

However, the results are far from this:

In[84]:= fit["BestFitParameters"]

Out[84]= {a -> 0.743405, b -> -0.203007}

sje1g13
  • 45
  • 4
  • 3
    I don't really understand your distinction between $x$ and $t$. You can rename one of them to the other one: ResourceFunction["MultiNonlinearModelFit"][Join[{datax, datat}], {modelx[a, b, x], modelt[a, b, x]}, {a, b}, {x}] Also, are you sure you actually want $a$ and $b$ to be the same in both models? Plotting your data, this clearly isn't the case ... – Domen Oct 09 '23 at 10:22
  • Thank you for your comments. You are certainly correct. I tried to simplify the problem at hand, but clearly went too far. I have edited my question to (hopefully) give a better idea of what I'd like to achieve. – sje1g13 Oct 09 '23 at 12:39
  • @sje1g13 You still try to fit using two different independent variables. You need to stick with just x or just t. Also: why are you adding the models together? – Sjoerd Smit Oct 09 '23 at 12:51
  • RandomReal[NormalDistribution[0, noiseLevel]] seems to be an odd construction. Might you want RandomVariate[NormalDistribution[0, noiseLevel]]. – JimB Oct 09 '23 at 16:07
  • The better approach is to use MultiNonlinearModelFit as @Domen suggests. But to make your code work try fit = NonlinearModelFit[Join[dataDec, dataPol], If[i == 1, modelDec[a, b, x], modelPol[a, b, x]], {{a, 1.0}, {b, 0.5}}, {i, x}, MaxIterations -> Infinity];. – JimB Oct 09 '23 at 17:01

1 Answers1

1

Your problem is that you want to fit 2 function simultaneously. Unfortunately, there is no fit function that can fit vector valued functions. Therefore, you need to "roll your own". You may achieve this by defining a function that calculates the sum of the squared errors and minimize this function.

Here is how you proceed:

The data:

nPoints = 256;
noiseLevel = 10^-1.5;

dataDecay = Table[{t, Sin[2 [Pi] 1.0 t] Exp[-0.5 t] + 0 RandomReal[NormalDistribution[0, noiseLevel]]}, {t, 0, 10, 10/nPoints}] // Chop;

dataPolys = Table[{x, 1 + (1.0 x^2 + 0.5 x)/100 + 0 RandomReal[NormalDistribution[0, noiseLevel]]}, {x, 0, 10, 10/nPoints}] // Chop;

The models:

modelDec[a_, b_, t_] := Sin[2 \[Pi] a t] Exp[-b t];
modelPol[a_, b_, x_] := 1 + (a x^2 + b x)/100;

The error function:

err[a_, b_] = 
  Total[(modelDec[a, b, #] & /@ dataDecay[[All, 1]] - 
       dataDecay[[All, 2]])^2 + (modelPol[a, b, #] & /@ 
        dataPolys[[All, 1]] - dataPolys[[All, 2]])^2];

And the minimization:

NMinimize[err[a, b], {a, b}]

{0.533559, {a -> 1.00172, b -> 0.503512}}

Daniel Huber
  • 51,463
  • 1
  • 23
  • 57