# Fitting different data sets with different models using FindFit

Posted 2 months ago
785 Views
|
10 Replies
|
2 Total Likes
|
 Referring to a fitting made previously thanks to the help of a member: Fit noisy data as a Least squares problem?I would like to ask for help to fit 2 noisy data sets to 2 models, referring to the same equation with different conditions. I have an equation T1: Dif = 0.00000013; K0 = 0.5; z = 0; linf\[Alpha] = 0.0; T1[x_, y_, t_, b_, l_, d_, m_, Intens_] /; NumberQ[x] && NumberQ[y] && NumberQ[b] && NumberQ[l] && NumberQ[d] && NumberQ[m] && NumberQ[Intens] := Intens/(2*\[Pi]*K0)* NIntegrate[ Sqrt[1 + m^2]* Erfc[Sqrt[(x - \[Alpha])^2 + (y - \[Beta])^2 + (z - (m*\[Alpha] - d))^2]/Sqrt[ Dif*4*t]]/(Sqrt[(x - \[Alpha])^2 + (y - \[Beta])^2 + (z - (m*\ \[Alpha] - d))^2]), {\[Beta], -b/2, b/2}, {\[Alpha], linf\[Alpha], l*(1/(1 + m^2))}] I create noisy data for y=0 and x=0: Block[{y = 0, t = 1, b = 0.001, l = 0.001, d = 0.0001, m = -Tan[45*Pi/180], Intens = 25000}, datax = Table[{x, T1[x, y, t, b, l, d, m, Intens] + Random[NormalDistribution[0, 0.1]]}, {x, -0.002, 0.002, 0.000015}]]; Block[{x = 0, t = 1, b = 0.001, l = 0.001, d = 0.0001, m = -Tan[45*Pi/180], Intens = 25000}, datay = Table[{y, T1[x, y, t, b, l, d, m, Intens] + Random[NormalDistribution[0, 0.1]]}, {y, -0.002, 0.002, 0.000015}]]; and after, all I want to do, is to fit the two noisy data sets to the equation above for y=0 and x=0, in order to find the best fitting parameters {b,l,d,m,Intens}.I thought about using Join as follow (t=1): alldata = Join[datax, datay]; fit=FindFit[alldata,T1[x,y,1,b,l,d,m,Intens],{b,l,d,m,Intens},{x,y},Method->"LevenbergMarquardt"] but nothing. Could anyone help me? Thanks.
10 Replies
Sort By:
Posted 2 months ago
 If I reformulate with Integrate instead of NIntegrate it does not give an error: FindFit[alldata, 0.3183098861837907 Intens Integrate[ Sqrt[1 + m^2]* Erfc[Sqrt[(x - \[Alpha])^2 + (y - \[Beta])^2 + (z - (m*\[Alpha] - d))^2]/Sqrt[ Dif*4*1]]/(Sqrt[(x - \[Alpha])^2 + (y - \[Beta])^2 + (z - (m*\ \[Alpha] - d))^2]), {\[Beta], -b/2, b/2}, {\[Alpha], linf\[Alpha], l*(1/(1 + m^2))}], {b, l, d, m, Intens}, {x, y}, Method -> "LevenbergMarquardt"] but it takes such a long time that I did not wait. The numerical evaluation of the double integral for so many different values of the parameters may be very onerous.
Posted 2 months ago
 I tried to perform this calculation, but I don't get any value from the parameters. I also get an error: FindFit::fitc: Number of coordinates (1) is not equal to the number of variables (2). 
Posted 2 months ago
 That is interesting.I had a look at your data ( datax and datay ) and tried this: ff[x_, a_, b_, c_] := a Exp[-b (x - c)^2] and sola = FindFit[datax, ff[x, a, b, c], {{a, 5.85}, {b, 5.1 10^6}, {c, 0.00015}}, x] (* {a -> 5.54423, b -> 4.52534*10^6, c -> 0.000158701} *) Then Plot[Evaluate[ff[x, a, b, c] /. sola], {x, -0.002, 0.002}, Epilog -> Point /@ datax, PlotRange -> {0, 7}] Similarly solb = FindFit[datay, ff[x, a, b, c], {{a, 5.05}, {b, 2.56 10^6}, {c, -0.00001}}, x] (* {a -> 5.54337, b -> 2.59886*10^6, c -> 1.4393*10^-7} *) Plot[Evaluate[ff[x, a, b, c] /. solb], {x, -0.002, 0.002}, Epilog -> Point /@ datay, PlotRange -> {0, 7}] These are the results from FindFit, and I believe it. Anyhow, one could get the impression (optically) that there is somethig better for datax, but note, this is an impression. Play around with a, b and c: Manipulate[ Plot[ff[x, a, b, c], {x, -0.002, 0.002}, PlotRange -> {0, 7}, Epilog -> Point /@ datax ], {a, 5, 7}, {b, .5 10^6., 10^7}, {c, -.002, .002}] 
Posted 2 months ago
 Thank you for the reply. Anyway, I don't understand well. ff in this case what does it refer to?
Posted 2 months ago
 Hello Lorenzo,ff is nothing but an empirical function to give an appoximation for your data.At least it shows that your complicated integral seems to have something to do with the error-function.And I think you can somewhat simplify your integral.If I didn't make a mistake your integrand is (could you please check this?) jj = (Sqrt[1 + m^2] Erfc[Sqrt[(x - \[Alpha])^2 + (d + z - m \[Alpha])^2 + (y - \[Beta])^2]/(2 Sqrt[Dif t])])/ Sqrt[(x - \[Alpha])^2 + (d + z - m \[Alpha])^2 + (y - \[Beta])^2]; Now introduce some abbreviations: rule = {k1 -> d^2 + x^2 + 2 d z + z^2 - (d m + x + m z)^2/(1 + m^2), k2 -> 1/(4 Dif t), \[Alpha]1 -> \[Alpha] - (d m + x + m z)/(1 + m^2), \[Beta]1 -> \[Beta] - y}; then you can write your integrand as jj3 = Sqrt[1 + m^2] Erfc[Sqrt[k1 k2 + k2 (1 + m^2) \[Alpha]1^2 + k2 \[Beta]1^2 ]]/ Sqrt[k1 + (1 + m^2) \[Alpha]1^2 + \[Beta]1^2]; jj - (jj3 /. rule) // Simplify ( I did not succeed to bring the arguments of the Erfc into the same form, but they are equal)I define a coordinate-transformation for the variables beta1 = fb and alpha1 = fa fb = Sqrt[(u + v)^2 - k1] Sin[v]; fa = Sqrt[((u + v)^2 - k1)/(1 + m^2)] Cos[v]; rulec = {\[Beta]1 -> fb, \[Alpha]1 -> fa}; and get as the new integrand jj4 = (jj3 /. rulec) Det[( { {D[fb, u], D[fb, v]}, {D[fa, u], D[fa, v]} } )] // PowerExpand // FullSimplify jj5 = jj4 /. Sqrt[k2 (u + v)^2] :> (u + v) Sqrt[k2] /. 1/Sqrt[(u + v)^2] :> 1/(u + v) namely -Erfc[Sqrt[k2] (u + v)] and this can be integrated explicitly Integrate[Integrate[-Erfc[Sqrt[k2] (u + v)], u], v] or Integrate[Integrate[-Erfc[Sqrt[k2] (u + v)], {u, u1, u2}], {v, v1, v2}] If you do it you will see that there are a couple of error - functions to be seen in the result.What I didn't do is to look for the boundaries for u and v to be compatible with the region spanned by alpha1 and beta1 (if this is possible at all).But I think you can do the backtransformation (with several solutions) sol = Solve[{fb == bb, fa == aa}, {u, v}]; and look for example at jj5 /. sol[[6]] /. rule or jj5 /. sol[[4]] /. rule 
Posted 2 months ago
 Now I get it. The problem is that I must necessarily use my original function. Thanks anyway.
Posted 2 months ago
 Summary: I want to fit the noisy data to the model I found them with (T1). The parameters are the same for both the fit to be made {b,l,d,m,Intens}. For datax the free variable is x. For datay the free variable is y. I don't want the problem to be due to this.
Posted 2 months ago
 You are fitting multivariate data (function dependent on x and y) but your data only has one coordinate logged (either x or y). so either add the respective x and y values to your data or just generate the full data range for x and y. Both are equally valid in this case. Dif = 0.00000013; K0 = 0.5; z = 0; linf\[Alpha] = 0.0; T1[x_, y_, t_, b_, l_, d_, m_, Intens_] /; NumberQ[x] && NumberQ[y] && NumberQ[b] && NumberQ[l] && NumberQ[d] && NumberQ[m] && NumberQ[Intens] := Intens/(2*\[Pi]*K0)* NIntegrate[ Sqrt[1 + m^2]* Erfc[Sqrt[(x - \[Alpha])^2 + (y - \[Beta])^2 + (z - (m*\[Alpha] - d))^2]/Sqrt[ Dif*4*t]]/(Sqrt[(x - \[Alpha])^2 + (y - \[Beta])^2 + (z - (m*\ \[Alpha] - d))^2]), {\[Beta], -b/2, b/2}, {\[Alpha], linf\[Alpha], l*(1/(1 + m^2))}] t = 1; b = 0.001; l = 0.001; d = 0.0001; m = -Tan[45*Pi/180]; Intens = 25000; step = 0.00006; datax = Table[{x, 0, T1[x, 0, t, b, l, d, m, Intens] + Random[NormalDistribution[0, 0.1]]}, {x, -0.002, 0.002, step}]; datay = Table[{0, y, T1[0, y, t, b, l, d, m, Intens] + Random[NormalDistribution[0, 0.1]]}, {y, -0.002, 0.002, step}]; alldata = Join[datax, datay]; dataxy = Flatten[ Table[{x, y, T1[x, y, t, b, l, d, m, Intens] + Random[NormalDistribution[0, 0.1]]}, {x, -0.002, 0.002, 4 step}, {y, -0.002, 0.002, 4 step}], 1]; Show[ Plot[{T1[x, 0, t, b, l, d, m, Intens], T1[0, x, t, b, l, d, m, Intens]}, {x, -0.002, 0.002}, PlotStyle -> Red], ListPlot[{datax[[All, {1, 3}]], datay[[All, {2, 3}]]}, PlotStyle -> {Black, Gray}] ] Show[ Plot3D[ T1[x, y, t, b, l, d, m, Intens], {x, -0.002, 0.002}, {y, -0.002, 0.002}, PlotRange -> Full, Mesh -> False, PlotStyle -> Gray, Lighting -> "Neutral"], ListPointPlot3D[dataxy, PlotRange -> Full, PlotStyle -> Black], Graphics3D[{Thick, Red, Line[datax], Line[datay]}] ] (*only lines*) fitL = FindFit[alldata, T1[x, y, 1, bf, lf, d, m, Intens], {{bf, 0.005}, {lf, 0.005}(*,d,m, Intens*)}, {x, y}, Method -> "LevenbergMarquardt"] (*full surface*) fitS = FindFit[dataxy, T1[x, y, 1, bf, lf, d, m, Intens], {{bf, 0.005}, {lf, 0.005}(*d,m, Intens*)}, {x, y}, Method -> "LevenbergMarquardt"] 
 I have another problem with the aforementioned fit. I would like to perform the fit considering also as an independent variable "t". The latter is described by two similar equations, the first (T1) one for 0 <= t <= tau and the second (T2) for t > tau. T2[x_, y_, t_, b_, l_, d_, tau_, m_, Intens_] /; NumberQ[x] && NumberQ[y] && NumberQ[t] && NumberQ[b] && NumberQ[l] && NumberQ[d] && NumberQ[tau] && NumberQ[m] && NumberQ[Intens] := Intens/(2*\[Pi]*K0)*Sqrt[1 + m^2]* NIntegrate[ Erfc[Sqrt[(x - \[Alpha])^2 + (y - \[Beta])^2 + (z - (m*\[Alpha] - d))^2]/Sqrt[ Dif*4*t]]/(Sqrt[(x - \[Alpha])^2 + (y - \[Beta])^2 + (z - (m*\ \[Alpha] - d))^2]) - Erfc[Sqrt[(x - \[Alpha])^2 + (y - \[Beta])^2 + (z - (m*\[Alpha] - d))^2]/Sqrt[ Dif*4*(t - tau)]]/(Sqrt[(x - \[Alpha])^2 + (y - \[Beta])^2 + (z - (m*\ \[Alpha] - d))^2]), {\[Beta], -b/2, b/2}, {\[Alpha], linf\[Alpha], l*(1/(1 + m^2))}] In producing noisy data I already get the first problems: Block[{x = 0, y = 0, b = 0.001, l = 0.001, d = 0.0001, m = -Tan[45*Pi/180], Intens = 25000}, datat1 = Table[{0, 0, t, T1[x, y, t, b, l, d, m, Intens] + Random[NormalDistribution[0, 0.1]]}, {t, 0, 1, 0.01}]]; NIntegrate::inumr: The integrand (Sqrt[2] Erfc[ComplexInfinity])/Sqrt[\[Alpha]^2+(0.0001 +\[Alpha])^2+\[Beta]^2] has evaluated to non-numerical values for all sampling points in the region with boundaries {{0,0.0005},{0.,0.0005}}. Power::infy: Infinite expression 1/0. encountered. General::stop: Further output of Power::infy will be suppressed during this calculation. Then I thought about doing the same thing for T2 and then using Join for the 4 data sets. But I think the problem starts already in the creation of noisy data.I also tried with: If[0 <= t <= 1, T1[0, 0, t, b, l, d, m, Intens], T2[0, 0, t, b, l, d, 1, m, Intens]] `but nothing happens.I wish I could fit the 4 datasets (for x, y and t) together with the 2 equations, always if this is possible.