Skip to main content

data - Fitting multiple datasets to the same nonlinear model


Perhaps a rather simple question, but the suggested questions seem to all be covering something slightly different. What I have can be described as follows.


One takes a dataset, call it Data; it has dimensions 10 by 1000 by 2 lets say, which is to say it is 10 arrays of {x,y} data of length 1000.



Now I know the following about it: it behaves according to a model Model[x,a,b,c]. Moreover, all 10 arrays have the same a,bbut not the same cparameter. What I am interested in is how one does a NonLinearModelFit over such a dataset, fitting a,b,c1...c10 in a single go. I should note that I have an initial value array cinit of dimensions 10 as well.


One way I thought of that should in principle work is with KroneckerDelta; I could prepend each dataset with some coordinate that I could use in addition to x and then make a model with the deltafunctions, but this feels very sloppy. There should be an easy way, should there not?


I would post the full model and the data I am working with, but not only is it a big dataset, the parameters and the model are still very much undetermined. I think it would make the question rather offtopic as it is my job to figure out the parameters, not you. So I thought it would make sense to just stick to the concept.



Answer



Here is a crude way to do it. First the dataset id is appended to the 1000 x 10 x 2 data array. Then the function you have is modified to create a dummy variable vector and the dot product is used to select the parameter of interest: c1, c2, c3,...., or c10.


(* Number of different values of c *)
nc = 10;
(* Number of observations per dataset *)
n = 1000; 


(* Create some data for a 1000 x 10 x 2 array *)
xx = Table[i/n, {i, n}];
data = Table[
   Transpose[{xx, 
     5 + 2 xx + xx^2 + RandomVariate[NormalDistribution[0, 1], n]}], {i, nc}];

(* Add in dataset number: 1 through nc *)
data2 = Flatten[
   Table[{i, data[[i, j, 1]], data[[i, j, 2]]}, {i, nc}, {j, n}], 1];


model[x_, a_, b_, i_, c_] := Module[{d, j},
   (* Create dummy variables *)
   d = Table[If[j == i, 1, 0], {j, Length[c]}];
   (* Determine function value *)
   a + b x + (c.d) x^2];

nlm = NonlinearModelFit[data2, 
   model[x, a, b, z, nc, {c1, c2, c3, c4, c5, c6, c7, c8, c9, c10}],
   {a, b, c1, c2, c3, c4, c5, c6, c7, c8, c9, c10}, {z, x}];
nlm["BestFitParameters"]

(* {a -> 4.99484, b -> 1.93603, c1 -> 1.18239, c2 -> 1.0829, 
 c3 -> 1.21282, c4 -> 1.06046, c5 -> 1.13936, c6 -> 1.07549, 
 c7 -> 1.01086, c8 -> 0.960816, c9 -> 0.979848, c10 -> 1.03909} *)

But note that I get a warning:



Experimental`NumericalFunction::dimsl: {x} given in {z,x} should be a list of dimensions for a particular argument.



Also, besides assuming common values for a and b , you are also assuming a common error variance. Residuals should be checked to see if there are any departures from that assumption or if a more complex error structure is warranted.


Comments

Post a Comment

Popular posts from this blog

mathematical optimization - Minimizing using indices, error: Part::pkspec1: The expression cannot be used as a part specification

I want to use Minimize where the variables to minimize are indices pointing into an array. Here a MWE that hopefully shows what my problem is. vars = u@# & /@ Range[3]; cons = Flatten@ { Table[(u[j] != #) & /@ vars[[j + 1 ;; -1]], {j, 1, 3 - 1}], 1 vec1 = {1, 2, 3}; vec2 = {1, 2, 3}; Minimize[{Total@((vec1[[#]] - vec2[[u[#]]])^2 & /@ Range[1, 3]), cons}, vars, Integers] The error I get: Part::pkspec1: The expression u[1] cannot be used as a part specification. >> Answer Ok, it seems that one can get around Mathematica trying to evaluate vec2[[u[1]]] too early by using the function Indexed[vec2,u[1]] . The working MWE would then look like the following: vars = u@# & /@ Range[3]; cons = Flatten@{ Table[(u[j] != #) & /@ vars[[j + 1 ;; -1]], {j, 1, 3 - 1}], 1 vec1 = {1, 2, 3}; vec2 = {1, 2, 3}; NMinimize[ {Total@((vec1[[#]] - Indexed[vec2, u[#]])^2 & /@ R...
Post a Comment
Read more

functions - Get leading series expansion term?

Given a function f[x] , I would like to have a function leadingSeries that returns just the leading term in the series around x=0 . For example: leadingSeries[(1/x + 2)/(4 + 1/x^2 + x)] x and leadingSeries[(1/x + 2 + (1 - 1/x^3)/4)/(4 + x)] -(1/(16 x^3)) Is there such a function in Mathematica? Or maybe one can implement it efficiently? EDIT I finally went with the following implementation, based on Carl Woll 's answer: lds[ex_,x_]:=( (ex/.x->(x+O[x]^2))/.SeriesData[U_,Z_,L_List,Mi_,Ma_,De_]:>SeriesData[U,Z,{L[[1]]},Mi,Mi+1,De]//Quiet//Normal) The advantage is, that this one also properly works with functions whose leading term is a constant: lds[Exp[x],x] 1 Answer Update 1 Updated to eliminate SeriesData and to not return additional terms Perhaps you could use: leadingSeries[expr_, x_] := Normal[expr /. x->(x+O[x]^2) /. a_List :> Take[a, 1]] Then for your examples: leadingSeries[(1/x + 2)/(4 + 1/x^2 + x), x] leadingSeries[Exp[x], x] leadingSeries[(1/x + 2 + (1 - 1/x...
Post a Comment
Read more

What is and isn't a valid variable specification for Manipulate?

I have an expression whose terms have arguments (representing subscripts), like this: myExpr = A[0] + V[1,T] I would like to put it inside a Manipulate to see its value as I move around the parameters. (The goal is eventually to plot it wrt one of the variables inside.) However, Mathematica complains when I set V[1,T] as a manipulated variable: Manipulate[Evaluate[myExpr], {A[0], 0, 1}, {V[1, T], 0, 1}] (*Manipulate::vsform: Manipulate argument {V[1,T],0,1} does not have the correct form for a variable specification. >> *) As a workaround, if I get rid of the symbol T inside the argument, it works fine: Manipulate[ Evaluate[myExpr /. T -> 15], {A[0], 0, 1}, {V[1, 15], 0, 1}] Why this behavior? Can anyone point me to the documentation that says what counts as a valid variable? And is there a way to get Manpiulate to accept an expression with a symbolic argument as a variable? Investigations I've done so far: I tried using variableQ from this answer , but it says V[1...
Post a Comment
Read more