2014-10-27 3 views
1

У меня есть функция, которая представляет собой произведение, составленное из трех (k) факторов. Каждый фактор представляет собой вычитание двух гауссовских CDF со случайными величинами R и L. Эти случайные величины определяются по 4 параметрам.Невозможно интегрировать вычитание двух гауссовских CDF [Wolfram Mathematica]

enter image description here

Код ниже показывает, как я построить основную функцию (в соответствии с двумя независимыми переменными д и е), и как случайные величины вычисляются

sigma = 1; 
k = 3; 
priors = {}; 
AppendTo[priors, 1/k + e]; 
Do[AppendTo[priors, 1/k - e/(k - 1)], {c, 2, k}]; 

L[priors_, sigma_, d_, i_] := Do[ 
maxVal = -Infinity; 
Do[ 
    val = (2*sigma^2*Log[priors[[i]]/priors[[j]]] + d^2 (j^2 - i^2 + 2 (i - j)))/(2 (j - i) d); 
    If[val > maxVal, maxVal = val, Null]; 
, {j, 1, i - 1}]; 
Return[maxVal]; 
, {1}]; 

R[priors_, sigma_, d_, i_] := Do[ 
minVal = Infinity; 
Do[ 
val = (2*sigma^2*Log[priors[[j]]/priors[[i]]] + d^2 (i^2 - j^2 + 2 (j - i)))/(2 (i - j) d); 
If[val < minVal, minVal = val, Null]; 
, {j, i + 1, k}]; 
Return[minVal]; 
, {1}]; 

Print[ 
Plot3D[ 
    Product[ 
    If[R[priors, sigma, d, c] < L[priors, sigma, d, c], 0, 
    (CDF[NormalDistribution[(c - 1) d, sigma], R[priors, sigma, d, c]] - 
    CDF[NormalDistribution[(c - 1) d, sigma], L[priors, sigma, d, c]])] 
    , {c, 1, k}] 
, {d, 0.01, 5} 
, {e, -1/k, 1 - 1/k}, PlotRange -> {All, All, All}, AxesLabel -> Automatic]]; 

Resulting plot of the first code

В настоящее время, Я хочу интегрировать функцию по d (в той же области, что и Plot3D, d = 0.01-5), и построить результаты по только независимой переменной e.

enter image description here

Ниже приведен код, который я использовал.

Print[ 
Plot[ 
    Integrate[ 
    Product[ 
    If[R[priors, sigma, d, c] < L[priors, sigma, d, c], 0, 
    (CDF[NormalDistribution[(c - 1) d, sigma], R[priors, sigma, d, c]] - 
     CDF[NormalDistribution[(c - 1) d, sigma], L[priors, sigma, d, c]])] 
    , {c, 1, k}] 
    , {d, 0.01, 5}] 
, {e, -1/k, 1 - 1/k}, PlotRange -> {All, All}, AxesLabel -> Automatic]]; 

Integrating over d

Однако полученный сюжет не то, что я ожидал. Это постоянно, и в 3D-графике видно, что этого не может быть. Кто-нибудь знает, что происходит и что делать, чтобы получить реальную интеграцию функции? Заранее спасибо.

ответ

2

При вычислении val внутри функций L и R результат является символическим (поскольку e не определен). Логический val < minVal Таким образом, неопределенный, и в результате minVal никогда не устанавливаются (так что L и R возвращение бесконечности каждый раз)

(очищено несколько других вещей ..)

sigma = 1; 
k = 3; 
priors = Join[ {1/k + e} , Table[1/k - e/(k - 1) , {c, 2, k} ] ]; 
L[priors0_, sigma_, d_, i_, e0_] := Module[{priors, maxVal, val, e}, 
    Do[maxVal = -Infinity; 
    priors = priors0 /. e -> e0 ; 
    Do[val = (2*sigma^2*Log[priors[[i]]/priors[[j]]] + 
     d^2 (j^2 - i^2 + 2 (i - j)))/(2 (j - i) d); 
     If[val > maxVal, maxVal = val];, {j, 1, i - 1}];, {1}]; maxVal]; 
R[priors0_, sigma_, d_, i_, e0_] := Module[{priors, maxVal, val, e}, 
    priors = priors0 /. e -> e0; 
    Min[Table[(2*sigma^2*Log[priors[[j]]/priors[[i]]] + 
     d^2 (i^2 - j^2 + 2 (j - i)))/(2 (i - j) d), {j, i + 1, k}]]]; 
g[d_?NumericQ, c_, e_] := 
    Product[If[R[priors, sigma, d, c, e] < L[priors, sigma, d, c, e], 
0, 
    (CDF[NormalDistribution[(c - 1) d, sigma], R[priors, sigma, d, c, e]] - 
    CDF[NormalDistribution[(c - 1) d, sigma], L[priors, sigma, d, c, e]])], 
    {c, 1, k}]; 
Plot[NIntegrate[g[d, c, e], {d, 0.01, 5}], {e, -1/k, 1 - 1/k}, 
     PlotRange -> {All, All}, AxesLabel -> Automatic] 

enter image description here

Смежные вопросы