0

最適化ライブラリ MPFIT を使用して、ガウス関数をデータに適合させようとしています。実際、このコードは MPFIT ライブラリに付属するサンプル コードの一部です。元のコードは、関数の導関数を数値的に内部的に自動的に計算し、完全に機能します。MPFIT ライブラリを使用すると、ユーザーは関数の派生物を提供することもできます。ここから問題が始まります。残差と関数の一次偏微分を計算するために使用される関数を次に示します。

int gaussfunc(int m, int n, double *p, double *dy, double **derivs, void *vars)
{
  int i,j;
  struct vars_struct *v = (struct vars_struct *) vars;
  double *x, *y, *ey;
  double a = p[1];
  double b = p[2];
  double c = p[3];
  double d = p[0];
  x = v->x;
  y = v->y;
  ey = v->ey;

  for (i=0; i<m; i++)
  {
    dy[i] = (y[i] - (a*exp(-(x[i]-b)*(x[i]-b)/(2*c*c))+d))/ey[i];
  }

  // the code below this point is the code I added to calculate the derivatives.
  if(derivs)
  {
      for(j = 0; j < n; j++)
      {
          if (derivs[j])
          {
              for (i = 0; i < m; i++)
              {
                  double da =     exp(-(x[i]-b)*(x[i]-b)/(2*c*c));
                  double db = a * exp(-(x[i]-b)*(x[i]-b)/(2*c*c)) * (x[i]-b)/(c*c);
                  double dc = a * exp(-(x[i]-b)*(x[i]-b)/(2*c*c)) * (x[i]-b)*(x[i]-b)/(c*c*c);
                  double dd = 1;
                  double foo;
                  if     (j == 0) foo = dd;
                  else if(j == 1) foo = da;
                  else if(j == 2) foo = db;
                  else if(j == 3) foo = dc;
                  derivs[j][i] = foo;
              }
          }
      }
  }

  return 0;
}

「if (derivs)」の行の上のコードは元のコードですが、リファクタリングされています。その下のコードは導関数を計算するための私のコードです。私の数学は正しいと信じており、https://math.stackexchange.com/questions/716545/calculating-the-first-order-partial-derivatives-of-the-gaussian-function/716553

ユーザー定義の導関数で MPFIT を使用しているときに同じ問題に遭遇した人はいますか?

ありがとうございました。

4

1 に答える 1

0

残差の計算は (DATA-MODEL)/SIGMA であるため、導関数は次のようになります。

[-d(モデル)/d(PARAM)]/シグマ。

したがって、この行：

derivs[j][i] = foo;

になります：

derivs[j][i] = -foo/ey[i];

問題が解決しました！ありがとう！

于 2014-03-19T02:03:30.483 に答える