1
$\begingroup$

My data gives some function integrated over bins of the argument. I want to fit the function with this data. I determine another function with the bin number as the argument and use NonlinearModelFit. But the bin number is passed to the function as a symbol so it fails:

Part::pkspec1: The expression jpp cannot be used as a part specification.

The (simplified) example follows. The data table contains lower and upper limits of integration for each bin, log10 of the measured value of the integral and its error. The function I really need is much more complicated, so the integral cannot be done analytically.

data = 
  {{-1.7211, -1.4201, 0.7799, 0.0076}, {-1.4201, -1.0221, 0.4844, 0.0105}, 
   {-1.0221, -0.545, -0.0118, 0.0184}, {-0.545, -0.0221, -0.9004, 0.0508}, 
   {-0.0221, 0.6768, -2.1792, 0.2902}};
y[a_, b_, jp_] := 
   Log[10., NIntegrate[a x + b, {x, data[[jp]][[1]], data[[jp]][[2]]}]];
nlm = 
  NonlinearModelFit[
    data[[All, 3]], y[a, b, jpp], {{a, 0.1}, b}, jpp, 
    Weights -> 1/data[[All, 4]]^2, VarianceEstimatorFunction -> (1 &)];

Using Extract instead of Part does not help.

I would appreciate your help.

Improve this question
$\endgroup$
3
  • $\begingroup$ That error is due to premature evaluation of your target function before jpp or jp have been assigned numerical values. You can prevent that by restricting y to be evaluated only with numerical arguments using NumericQ, I.e. changing its definition to Clear[y]; y[a_?NumericQ, b_?NumericQ, jp_?NumericQ] := yourcode. See this FAQ. $\endgroup$
    – MarcoB
    Commented Jun 27, 2018 at 13:00
  • $\begingroup$ Thank you MarcoB. This removes the error message. $\endgroup$
    – otfinr
    Commented Jun 27, 2018 at 13:29
  • $\begingroup$ There are things to do after your question is answered. It's a good idea to stay vigilant for some time, better approaches may come later improving over previous replies. Experienced users may point alternatives, caveats or limitations. New users should test answers before voting and wait 24 hours before accepting the best one. Participation is essential for the site, please do your part. $\endgroup$
    – rhermans
    Commented Jul 4, 2018 at 7:41

2 Answers 2

Reset to default
1
$\begingroup$

Your definition of y requires its 3rd argument to be a positive integer (even though you don't specify that constraint), but NonlinearModelFit does not support such models, which can be seen from 

data =
  {{-1.7211, -1.4201, 0.7799, 0.0076}, {-1.4201, -1.0221, 0.4844, 
    0.0105}, {-1.0221, -0.545, -0.0118, 
    0.0184}, {-0.545, -0.0221, -0.9004, 0.0508}, {-0.0221, 
    0.6768, -2.1792, 0.2902}};

Clear[y]
y[a_?NumericQ, b_?NumericQ, jp_Integer?Positive] := 
  Log[10, NIntegrate[a x + b, {x, data[[jp]][[1]], data[[jp]][[2]]}]]

nlm = 
  NonlinearModelFit[
    data[[All, 3]], {y[a, b, jpp], {jpp ∈ Integers}}, {{a, 0.1}, b}, jpp, 
    Weights -> 1/data[[All, 4]]^2, 
    VarianceEstimatorFunction -> (1 &)

msg

Improve this answer
$\endgroup$
1
  • $\begingroup$ Thank you m_goldberg. Requiring jp to be numeric (not necessarily integer), as suggested by MarcoB, solved the problem (at least it works without error messages). Hope I can trust the fit result (will check). Is there any other way to do the fit I need with Mathematica, without writing the full chi-squared routine? $\endgroup$
    – otfinr
    Commented Jun 27, 2018 at 13:35
0
$\begingroup$

With the actual function (that you are apparently not able to share) and/or the actual data, you'll either need to add in appropriate constraints for the parameters (as the value in the Log function will need to be positive) or have awfully good starting values.

Below is one approach to obtain the constraints and rewrite the function so that you won't have to worry about the issue that @m_goldberg pointed out:

constraints = Reduce[# > 0 & /@ (Integrate[a x + b, {x, #[[1]], #[[2]]}] & /@ data)]
(* b > 0 && -3.05483 b < a < 0.636699 b *)

y[a_?NumericQ, b_?NumericQ, x1_?NumericQ, x2_?NumericQ] := 
  Log[10., NIntegrate[a x + b, {x, x1, x2}]];
nlm = NonlinearModelFit[data[[All, {1, 2, 3}]], {y[a, b, x1, x2], constraints},
  {{a, 0.1}, b}, {x1, x2}, Weights -> 1/data[[All, 4]]^2, 
  VarianceEstimatorFunction -> (1 &)];
nlm["BestFitParameters"]
(* {a -> -6.73818, b -> 2.21297} *)
Improve this answer
$\endgroup$
1
  • $\begingroup$ Thank you JimB, it works for the toy example. It even worked for the actual function (which is an InterpolatingFunction) but took eternity to Reduce. I can probably try to have sufficiently good starting values. $\endgroup$
    – otfinr
    Commented Jun 29, 2018 at 14:19

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.