Tôi muốn tích hợp hàm mật độ xác suất từ (-\infty, a]
do không có sẵn tệp cdf ở dạng đóng. Nhưng tôi không chắc chắn làm thế nào để làm điều này trong C + +.Thường trình bậc hai cho mật độ xác suất
Nhiệm vụ này khá đơn giản trong Mathematica; Tất cả những gì tôi cần làm là xác định hàm,
f[x_, lambda_, alpha_, beta_, mu_] :=
Module[{gamma},
gamma = Sqrt[alpha^2 - beta^2];
(gamma^(2*lambda)/((2*alpha)^(lambda - 1/2)*Sqrt[Pi]*Gamma[lambda]))*
Abs[x - mu]^(lambda - 1/2)*
BesselK[lambda - 1/2, alpha Abs[x - mu]] E^(beta (x - mu))
];
và sau đó gọi số NIntegrate
Định kỳ để tích hợp số.
F[x_, lambda_, alpha_, beta_, mu_] :=
NIntegrate[f[t, lambda, alpha, beta, mu], {t, -\[Infinity], x}]
Bây giờ tôi muốn đạt được điều tương tự trong C++. Tôi sử dụng thường lệ gsl_integration_qagil
từ thư viện số gsl. Nó được thiết kế để tích hợp các chức năng trong khoảng thời gian bán vô hạn (-\infty, a]
mà chỉ là những gì tôi muốn. Nhưng tiếc là tôi không thể làm cho nó hoạt động được.
Đây là hàm mật độ trong C++,
density(double x)
{
using namespace boost::math;
if(x == _mu)
return std::numeric_limits<double>::infinity();
return pow(_gamma, 2*_lambda)/(pow(2*_alpha, _lambda-0.5)*sqrt(_pi)*tgamma(_lambda))* pow(abs(x-_mu), _lambda - 0.5) * cyl_bessel_k(_lambda-0.5, _alpha*abs(x - _mu)) * exp(_beta*(x - _mu));
}
Sau đó, tôi cố gắng và tích hợp để có được những lũy bằng cách gọi thông thường GSL.
cdf(double x)
{
gsl_integration_workspace * w = gsl_integration_workspace_alloc (1000);
double result, error;
gsl_function F;
F.function = &density;
double epsabs = 0;
double epsrel = 1e-12;
gsl_integration_qagil (&F, x, epsabs, epsrel, 1000, w, &result, &error);
printf("result = % .18f\n", result);
printf ("estimated error = % .18f\n", error);
printf ("intervals = %d\n", w->size);
gsl_integration_workspace_free (w);
return result;
}
Tuy nhiên gsl_integration_qagil
trả về lỗi, number of iterations was insufficient
.
double mu = 0.0f;
double lambda = 3.0f;
double alpha = 265.0f;
double beta = -5.0f;
cout << cdf(0.01) << endl;
Nếu tôi tăng kích thước của không gian làm việc thì chức năng bessel sẽ không đánh giá.
Tôi đã tự hỏi liệu có ai có thể cho tôi bất kỳ thông tin chi tiết nào về vấn đề của tôi hay không. Một cuộc gọi đến hàm Mathematica tương ứng F ở trên với x = 0.01
trả về 0.904384
.
Có thể mật độ được tập trung xung quanh một khoảng thời gian rất nhỏ (ví dụ: bên ngoài [-0.05, 0.05]
mật độ gần như là 0
, một ô được đưa ra dưới đây). Nếu có thể làm gì về điều này. Cảm ơn vì đã đọc.
Chức năng này có vẻ đối xứng, nghĩa là 'cdf (0) = 1/2'.Hãy nhớ rằng cdf được đánh giá tại x là giống như tích phân từ 0 đến x cộng với cdf được đánh giá ở 0. Tất nhiên, tôi đang đi từ hình dạng của biểu đồ, nó có thể không thực sự chính xác đối xứng. –
Xin chào @Ben, tôi ước rằng đó là trường hợp! – mark