Message Boards Message Boards

0
|
7312 Views
|
11 Replies
|
4 Total Likes
View groups...
Share
Share this post:

Find Fit a function

Posted 5 years ago

Hello, I'm trying to use Mathematica to fit a function, but have an error, can someone tell me what the problem is?

I first define the data table that will be fitted.

xmdatapar = 
  ReadList["/home/roberto/F?ica/Pesquisa/<X^m>/kaon/lib/pares/xm_\
pares.dat", {Number}];
xmpardata = Table[{xmdatapar[[i, 1]]}, {i, 1, Length[xmdatapar]}]

{{0.243098}, {0.127303}, {0.0834994}, {0.0609518}, {0.0473757}, \
{0.0384113}, {0.0320716}, {0.027407}, {0.0237647}, {0.020894}, \
{0.0185901}, {0.0167308}, {0.0151393}, {0.013839}, {0.0127182}, \
{0.0117066}, {0.0108577}, {0.0101237}, {0.00948312}, {0.00887854}, \
{0.00834716}, {0.00788257}, {0.00746653}, {0.00705183}, {0.00669624}}

in the second place the function with the variables that I want to fit.

\[Phi]par[x_] := 
  Npar*(x (1 - x))^
   aRL\[Alpha]par (1 + \[Alpha]2RLpar*
      GegenbauerC[2, \[Alpha]RLpar, (2 x - 1)]);

and finally I make the fit

fitpar = FindFit[
  xmpardata, \[Phi]par, {Npar, aRL\[Alpha]par, \[Alpha]2RLpar}, x]

FindFit::fitc: Number of coordinates (0) is not equal to the number of variables (1).

FindFit[{{0.243098}, {0.127303}, {0.0834994}, {0.0609518}, \
{0.0473757}, {0.0384113}, {0.0320716}, {0.027407}, {0.0237647}, \
{0.020894}, {0.0185901}, {0.0167308}, {0.0151393}, {0.013839}, \
{0.0127182}, {0.0117066}, {0.0108577}, {0.0101237}, {0.00948312}, \
{0.00887854}, {0.00834716}, {0.00788257}, {0.00746653}, {0.00705183}, \
{0.00669624}}, \[Phi]par, {Npar, aRL\[Alpha]par, \[Alpha]2RLpar}, x]

Thanks.

Attachments:
POSTED BY: RCS LV
11 Replies
Posted 5 years ago

ohh i see , but why normf[n_, a_, a2_] := Quiet@Norm@(xd[n, a, a2] - xmpardata)? thanks!

POSTED BY: RCS LV

Norm is a way to estimate the sum of the residuals and Quiet stops the beeping that I was getting from warnings. Here is an example:

dat = {d1, d2, d3};
model = {m1, m2, m3};
resid = model - dat;
Norm@resid // InputForm (* \
Sqrt[Abs[-d1+m1]^2+Abs[-d2+m2]^2+Abs[-d3+m3]^2] *)
POSTED BY: Tim Laska

Roberto,

I was able to get FindFit to work by using LocalAdaptive stepping on the NIntegrate and by setting the WorkingPrecision to 8.

data = Table[{2 i, xmpardata[[i]]}, {i, 1, 25}];
phi[n_, a_, a2_][x_] := 
 n*(x (1 - x))^(a - 1/2) (1 + a2*GegenbauerC[2, a, (2 x - 1)])
fitxm[n_, a_, a2_][m_] := 
 Quiet@NIntegrate[phi[n, a, a2][x] (2*x - 1)^m, {x, 0, 1}, 
   Method -> "LocalAdaptive"]
xd[fit_] := {#, (fitxm[n, a, a2] /. fit)[#]} & /@ Range[2, 50, 2]
title[fit_] := 
 StringTemplate["n = `n`\n a = `a`\n a2 = `a2`"][
  AssociationThread @@ ({(ToString@#) & /@ Keys@fit, Values@fit})]
fit = FindFit[data, 
   fitxm[n, a, a2][m], {{n, 1.1765`}, {a, 0.58}, {a2, -0.44}}, m, 
   WorkingPrecision -> 8];
ListPlot[{data, xd[fit]}, PlotLabel -> title[fit], 
 PlotLegends -> {"Data", "Fit"}]

Fitted

POSTED BY: Tim Laska
Posted 5 years ago

Thanks a lot for helping me!

I fix the code, but the fit is not good.

1. I read the data here my x

xmdatapar = 
  ReadList["/home/roberto/F?ica/Pesquisa/<X^m>/kaon/lib/pares/xm_par.\
dat", {Number}];
xmpardata = Table[xmdatapar[[i, 1]], {i, 1, Length[xmdatapar]}]
{0.243098, 0.127303, 0.0834994, 0.0609518, 0.0473757, 0.0384113, \
0.0320716, 0.027407, 0.0237647, 0.020894, 0.0185901, 0.0167308, \
0.0151393, 0.013839, 0.0127182, 0.0117066, 0.0108577, 0.0101237, \
0.00948312, 0.00887854, 0.00834716, 0.00788257, 0.00746653, \
0.00705183, 0.00669624}

2. i need to fit this variables *Npar, (\[Alpha]RLpar , \[Alpha]2*

\[Phi]par[x_] := 
 Npar*(x (1 - x))^(\[Alpha]RLpar - 
    1/2 ) (1 + \[Alpha]2*GegenbauerC[2, \[Alpha]RLpar, (2 x - 1)]);

3.I did a model to fit

Table[{m, Integrate[\[Phi]par[x] ( 2*x - 1)^m, {x, 0, 1}]}, {m, 2, 6, 
   2}];

fitxm = Integrate[\[Phi]par[x] ( 2*x - 1)^m, {x, 0, 1}];
  1. where m is numbers even 2,4,6...50 and x my moments, [Alpha]2 need to be a negative low value like around -0.30~-0.50

      FindFit[xmpardata, fitxm, {\[Alpha]RLpar, \[Alpha]2, Npar}, m]
    
        FindFit::nrlnum: The function value {-0.243098+0. I,-0.000149731-1.55718*10^-17 I,-0.0834994+0. I,0.00219318 -1.54661*10^-17 I,-0.0473757+0. I,0.000892465 -1.444*10^-17 I,-0.0320716+0. I,<<11>>,-0.00948312+0. I,-0.000897592-9.77385*10^-18 I,-0.00834716+0. I,-0.00090295-3.41989*10^-17 I,-0.00746655+0. I,-0.000880321-9.0694*10^-18 I,-0.00669639+0. I} is not a list of real numbers with dimensions {25} at {\[Alpha]RLpar,\[Alpha]2,Nc} = {1.00001,-0.0133867,1.31278}.
    
        {\[Alpha]RLpar -> 1.00001, \[Alpha]2 -> -0.0133867, Npar -> 1.31278}
        ListPlot[xmpardata]
        In[8]:= Table[
         fitxm /. {\[Alpha]RLpar -> 
            1.00001204228849, \[Alpha]2 -> -0.0133867063850715, 
           Npar -> 1.3127808940632832}, {m, 2, 50, 2}]
    
        Out[8]= {0.127153, 0.063145, 0.0393038, 0.0274371, 0.0205373, \
        0.0161126, 0.0130764, 0.010887, 0.00924693, 0.00798095, 0.00697962, \
        0.00617151, 0.00550909, 0.00495871, 0.00444164, 0.00412147, \
        -0.000131836, -0.036299, 0.1193, -0.170002, -4.57313, -135.358, \
        1499.52, 4139.2, 44729.1}
    
        In[10]:= ListPlot[
         Table[fitxm /. {\[Alpha]RLpar -> 
             1.00001204228849, \[Alpha]2 -> -0.0133867063850715, 
            Npar -> 1.3127808940632832}, {m, 2, 50, 2}]]
    
POSTED BY: RCS LV

Roberto,

I think that your integral expression leads to difficult numerics. If you take the semicolon away, the output is conditional expression about a page long. If you evaluate numerically, you often can get warnings and imaginary components in the results.

You may want use Mathematica's Manipulate functionality to manually explore the data. Here is an example where you try to manually minimize the Norm function by manipulating the sliders as shown below:

phi[n_, a_, a2_][x_] := 
 n*(x (1 - x))^(a - 1/2) (1 + a2*GegenbauerC[2, a, (2 x - 1)])
fitxm[n_, a_, a2_][m_] := 
 Quiet@NIntegrate[phi[n, a, a2][x] (2*x - 1)^m, {x, 0, 1}]
xd[n_, a_, a2_] := fitxm[n, a, a2][#] & /@ Range[2, 50, 2]
normf[n_, a_, a2_] := Quiet@Norm@(xd[n, a, a2] - xmpardata)
Manipulate[global = {n, a, a2}; 
 ListPlot[{xmpardata, xd[n, a, a2]}, PlotLabel -> normf[n, a, a2], 
  PlotLegends -> {"Data", "Fit"}], {{n, 1.1}, 0.9, 1.2}, {{a, 0.58}, 
  0.45, 0.7}, {{a2, -0.44}, -0.5, -0.3}]
Dynamic@global

Manipulate Example

POSTED BY: Tim Laska

If you look closely at you expression for $\phi par$, then you will see that you have 4 parameters versus the 3 specified in FindFit. One of them is most likely a typo. You can check the unique symbols in your expression as shown below:

\[Phi]par := 
  Npar*(x (1 - x))^
    aRL\[Alpha]par (1 + \[Alpha]2RLpar*
      GegenbauerC[2, \[Alpha]RLpar, (2 x - 1)]);
DeleteDuplicates@Cases[\[Phi]par, _Symbol, Infinity]
(* {Npar, x, aRL\[Alpha]par, \[Alpha]2RLpar, \[Alpha]RLpar} *)

enter image description here

POSTED BY: Tim Laska
Posted 5 years ago

oh i see , the right is

\[Alpha]RLpar = `aRL\[Alpha]par-1/2
POSTED BY: RCS LV

If your end goal is to perform the numerical integration, then you are probably making it much more difficult than it needs to be. If you plot the data with ListLogLogPlot, it is very linear suggesting LinearModelFit is in order.

data = {{2, 0.243098}, {4, 0.127303}, {6, 0.0834994}, {8, 
    0.0609518}, {10, 0.0473757}, {12, 0.0384113}, {14, 
    0.0320716}, {16, 0.027407}, {18, 0.0237647}, {20, 0.020894}, {22, 
    0.0185901}, {24, 0.0167308}, {26, 0.0151393}, {28, 0.013839}, {30,
     0.0127182}, {32, 0.0117066}, {34, 0.0108577}, {36, 
    0.0101237}, {38, 0.00948312}, {40, 0.00887854}, {42, 
    0.00834716}, {44, 0.00788257}, {46, 0.00746653}, {48, 
    0.00705183}, {50, 0.00669624}};
fit = LinearModelFit[Log@data, x, x]
phipar[x_] := Exp@Normal@fit;
fit["ParameterTable"]
Plot[Normal@fit, {x, 0.5, 4}, Axes -> False, Frame -> True, 
 Epilog -> {PointSize@.02, Point@Log@data}]
f[m_] := NIntegrate[(2*x - 1)^m (phipar[x]), {x, 0, 1}]
f[2] (* 0.130834 *)

Simple Approach

With your current formulation, $aRL\alpha par$ must be even integers to avoid complex numbers since $x\geq 2$ making the argument negative.

$$\phi_{par}=N_{par}\left ( x\left ( 1-x \right ) \right )^{aRL\alpha par}....$$

POSTED BY: Tim Laska
Posted 5 years ago

Hi Roberto,

There are two problems in your code:

  1. The table xmpardata needs to contain pairs which describe the value of x and the value of the model evaluated at x. Your table contains lists of only one element.
  2. The FindFit function wants a model which is a symbolic form. You are passing it a function phiPar. If you need the function, you can always construct the model by evaluating phiPar[x], otherwise, just define the model as the expression which defines phiPar.

Kind regards,

David

POSTED BY: David Keith
Posted 5 years ago

Sorry, I've never used Mathematica before, so I do not know many things, I have a dat with the first 25 even Gegenbauer polynomials of second order for fit, I fix the table

xmdatapar = 
  ReadList["/home/roberto/F?ica/Pesquisa/<X^m>/kaon/lib/pares/xm_\
pares.dat", {Number, Number}];
xmpardata = 
 Table[{xmdatapar[[i, 1]], xmdatapar[[i, 2]]}, {i, 1, 
   Length[xmdatapar]}]

{{2, 0.243098}, {4, 0.127303}, {6, 0.0834994}, {8, 0.0609518}, {10, 
  0.0473757}, {12, 0.0384113}, {14, 0.0320716}, {16, 0.027407}, {18, 
  0.0237647}, {20, 0.020894}, {22, 0.0185901}, {24, 0.0167308}, {26, 
  0.0151393}, {28, 0.013839}, {30, 0.0127182}, {32, 0.0117066}, {34, 
  0.0108577}, {36, 0.0101237}, {38, 0.00948312}, {40, 
  0.00887854}, {42, 0.00834716}, {44, 0.00788257}, {46, 
  0.00746653}, {48, 0.00705183}, {50, 0.00669624}}

where the first value is the first term sum of even numbers 2,4,6,8,10...50 and the second value is a x moments,

i need to fit this function

 In[86]:= \[Phi]par := 
      Npar*(x (1 - x))^
       aRL\[Alpha]par (1 + \[Alpha]2RLpar*
          GegenbauerC[2, \[Alpha]RLpar, (2 x - 1)]);

where phi is related to this integral

F[m_] := NIntegrate[(2*x - 1)^m (\[Phi]par[x] ), {x, 
   0, 1}]

where m = 2,4,6,8,10,12... so how I can do this ??

In[89]:= fitpar = 
 FindFit[xmpardata, {\[Phi]par}, {Npar, 
   aRL\[Alpha]par, \[Alpha]2RLpar}, {x}]

During evaluation of In[89]:= FindFit::nrlnum: The function value {-0.243098-2. (1. +1. (-1. \[Alpha]RLpar+18. \[Alpha]RLpar Plus[<<2>>])),-0.127303-12. (1. +1. (-1. \[Alpha]RLpar+98. \[Alpha]RLpar Plus[<<2>>])),<<22>>,-0.00669624-2450. (1. +1. (-1. \[Alpha]RLpar+19602. \[Alpha]RLpar Plus[<<2>>]))} is not a list of real numbers with dimensions {25} at {Npar,aRL\[Alpha]par,\[Alpha]2RLpar} = {1.,1.,1.}.

During evaluation of In[89]:= FindFit::nrlnum: The function value {-0.243098-2. (1. +1. (-1. \[Alpha]RLpar+18. \[Alpha]RLpar Plus[<<2>>])),-0.127303-12. (1. +1. (-1. \[Alpha]RLpar+98. \[Alpha]RLpar Plus[<<2>>])),<<22>>,-0.00669624-2450. (1. +1. (-1. \[Alpha]RLpar+19602. \[Alpha]RLpar Plus[<<2>>]))} is not a list of real numbers with dimensions {25} at {Npar,aRL\[Alpha]par,\[Alpha]2RLpar} = {1.,1.,1.}.

Out[89]= FindFit[{{2, 0.243098}, {4, 0.127303}, {6, 0.0834994}, {8, 
   0.0609518}, {10, 0.0473757}, {12, 0.0384113}, {14, 0.0320716}, {16,
    0.027407}, {18, 0.0237647}, {20, 0.020894}, {22, 0.0185901}, {24, 
   0.0167308}, {26, 0.0151393}, {28, 0.013839}, {30, 0.0127182}, {32, 
   0.0117066}, {34, 0.0108577}, {36, 0.0101237}, {38, 
   0.00948312}, {40, 0.00887854}, {42, 0.00834716}, {44, 
   0.00788257}, {46, 0.00746653}, {48, 0.00705183}, {50, 
   0.00669624}}, {Npar ((1 - x) x)^
   aRL\[Alpha]par (1 + \[Alpha]2RLpar (-\[Alpha]RLpar + 
        2 (-1 + 2 x)^2 \[Alpha]RLpar (1 + \[Alpha]RLpar)))}, {Npar, 
  aRL\[Alpha]par, \[Alpha]2RLpar}, {x}]

Thanks.

POSTED BY: RCS LV
Posted 5 years ago

While I think that datasets should explicitly have the response value and all predictors, it seems that most if not all of the Wolfram Language fitting functions take a dataset with just the response values in a single list and then assumes that the single predictor variable is from the list Range[Length[data]].

However, the OP's data would need to be Flatten'ed to get it in the form of a single list that those functions would accept.

Also having the term in the form of (x(1-x))^k and all positive response values suggests that maybe the x values are between zero and one rather than between 2 and 50.

POSTED BY: Jim Baldwin
Reply to this discussion
Community posts can be styled and formatted using the Markdown syntax.
Reply Preview
Attachments
Remove
or Discard

Group Abstract Group Abstract