MATLAB 이상한 오류 감마 함수 수치 적분

Question

i는 수치 적으로 통합하기 위해 다음과 같은 실행하려고 :

TPDF는 MATLAB 기능과는 skewtdis입니다

nu = 8; 
psi=-0.2; 
lambda = 1; 
git = @(u) tpdf((0 - lambda * skewtdis_inverse(u, nu, psi)), nu); 
g(t,i) = integral(git,1e-10,1-1e-10,'AbsTol',1e-16); 

:

function inv = skewtdis_inverse(u, nu, lambda) 
% PURPOSE: returns the inverse cdf at u of Hansen's (1994) 'skewed t' distribution 

c = gamma((nu+1)/2)/(sqrt(pi*(nu-2))*gamma(nu/2)); 
a = 4*lambda*c*((nu-2)/(nu-1)); 
b = sqrt(1 + 3*lambda^2 - a^2); 

if (u<(1-lambda)/2); 
inv = (1-lambda)/b*sqrt((nu-2)./nu)*tinv(u/(1-lambda),nu)-a/b; 
elseif (u>=(1-lambda)/2); 
inv = (1+lambda)/b*sqrt((nu-2)./nu).*tinv(0.5+1/(1+lambda)*(u-(1-lambda)/2),nu)-a/b; 
end 

: 역은 다음과 같습니다 내가 꺼내는 것은 다음과 같습니다 :

c = 감마 ((nu + 1)/2)/(sqrt (pi * (nu-2)) * 감마 (nu/2) ;

"F : \ Xyz \ skewtdis_inverse.m> skewtdis_inverse"를 호출하는 동안 출력 인수 "inv"(및 기타 항목)가 지정되지 않았습니다. @에서

오류 (U) TPDF ((0, λ * skewtdis_inverse (U, 민, PSI)), 뉴) integralCalc/iterateScalarValued (라인 314) FX = FUN (t)에

오류;

integralCalc/vadapt (line 133)의 오류 [q, errbnd] = iterateScalarValued (u, tinterval, pathlen);

integralCalc (line 76)의 오류 [q, errbnd] = vadapt (@ AtoBInvTransform, interval);

정수 중 오류 (줄 89) Q = integralCalc (fun, a, b, opstruct); (- 람다 * skewtdis_inverse (1E-10 민, PSI)), ν (0)

TPDF : 나 그러나 THR의 기능 처리를 호출하면

직접 아무런 문제가 없다

ANS =

1.4092e-11

01 23,516,

TPDF ((0 - 람다 * skewtdis_inverse (1-1e-10 민, PSI)), 뉴)

의 ANS =

7.0108e-10

내 노력은 매우 높이 평가됩니다!

Answer 1

기본적으로 integral은 함수 핸들이 벡터 입력을 필요로합니다. u의 모든 요소가 인 경우에만 조건은 true으로 계산되므로 코드에서 if- 구문은 복잡성을 생성합니다. u이 (1-lambda)/2보다 크고 작은 요소를 포함하는 벡터 인 경우 inv은 할당되지 않습니다.

1. 가 for -loop에 if -statement을 넣고 u의 모든 요소를 ​​반복 :

   두 가지 옵션이 있습니다.

2. 할당에 논리 인덱스를 사용합니다.

두 번째 옵션은 내 의견, 청소기, 빨리 큰 요소 수에 대한입니다 :

inv = u; % Allocation 

IsBelow = u < (1-lambda)/2; % Below the threshold 
IsAbove = ~IsBelow   ; % Above the threshold 

inv(IsBelow) = (1-lambda)/b*sqrt((nu-2)./nu)*tinv(u(IsBelow)/(1-lambda),nu)-a/b; 
inv(IsAbove) = (1+lambda)/b*sqrt((nu-2)./nu)*tinv(0.5+1/(1+lambda)*(u(IsAbove)-(1-lambda)/2),nu)-a/b;