Problem in transition matrix #189

khoilx · 2024-10-10T23:38:47Z

Hello Prof. Maih,
I'm facing the problem of retcode "3" (Problem in transition matrix) because I declared an alien list to calculate the transitional probability (tp_1_2) from an outside function (that calculate the mvncdf of a vector) (e.g. ! s1_tp_1_2 =tp_func(...)). Could you please help me to deal with this problem in RISE?

jmaih · 2024-10-11T12:22:58Z

Hi Khoi, Are you able to evaluate the function in matlab? Do you pass the alien function correctly? e.g. m=set(m, 'alien_list',{'yourFunctionName'}); J.

…

On Fri, Oct 11, 2024 at 1:39 AM khoilx ***@***.***> wrote: Hello Prof. Maih, I'm facing the problem of retcode "3" (Problem in transition matrix) because I declared an alien list to calculate the transitional probability (tp_1_2) from an outside function (that calculate the mvncdf of a vector) (e.g. ! s1_tp_1_2 =tp_func(...)). Could you please help me to deal with this problem in RISE? — Reply to this email directly, view it on GitHub <#189>, or unsubscribe <https://github.com/notifications/unsubscribe-auth/AATKBT3XWSHIUBCXSDIFB6DZ24FZXAVCNFSM6AAAAABPX2BP4GVHI2DSMVQWIX3LMV43ASLTON2WKOZSGU4DAMBWGI3DMNQ> . You are receiving this because you are subscribed to this thread.Message ID: ***@***.***>

khoilx · 2024-10-11T14:51:58Z

Hello Prof. Maih,
Yes, I pass the alien function correctly (I declared as:
m=rise('cmr_mod','steady_state_file','model_ss',...
'rise_flags',model_settings,...
'alien_list','tp_func2');).
I'm sending you the 'tp_func2.m' function file. Could you please check it for me?

function y=tp_func2(z1,z2,z3,a1,a2,a3,xmin,xmax,ymin,ymax,zmin,zmax,wmin,wmax,varargin)

if ~isreal(z1) && ~isreal(z2) && ~isreal(z3)

iz1=imag(z1);

iz2 = imag(z2);

iz3 = imag(z3);

z1=real(z1);

z2=real(z2);

z3=real(z3);

if abs(iz1)>1e-12 && abs(iz2)>1e-12 && abs(iz3)>1e-12 
    
    y=nan;
    
    return
    
end

end

Pev_sq = [a1 a2; a2 a3];

Pvv = [1, z3; z3, 1];

A = [z1 0; 0 z2];

try
% Try solving Sigw = ASigwA'+Pvv for Sigw.
Sigw = dlyap(A,Pvv);
catch
% If A is not stable, throw error and skip all calculations.
warning('A is not stable. Assigning 0 to likelihood');
Sigw = eye(2)*.1e10; %replace s.Sigw by garbage
%lik = 0; %replace s.lik by 0
end

%compute the covariance Pvve
Pvve = Pvv - Pev_sq;

%compute covariance matrix for 4d multivariate normal in the calculation of time-varying transition
Sigma_mvn = [Sigw, SigwA'; ASigw, Pvve + ASigwA'];

if isempty(varargin)

y=main();

else

if length(varargin)==2
    
    if ~strcmp(varargin{1},'diff')
        
        error('Wrong call to function')
        
    end
    
    switch varargin{2}
        
        case 1
            
            y=linhbin_z();
            
        case {2,3,4,5}
            
            error('Higher-order derivatives not implemented')
            
    end
    
else
    
    error('Wrong call to function')
    
end

end

function y=main()
    
    y=mvncdf([xmin,ymin,zmin,wmin],[xmax,ymax,zmax,wmax],[0 0 0 0],Sigma_mvn)/mvncdf([zmin,wmin],[zmax,wmax],[0 0],Sigw);
    
end

function fz=linhbin_z()
    
    fz=(mvnpdf([xmax,ymax,zmax,wmax],[0 0 0 0],Sigma_mvn)-mvnpdf([xmin,ymin,zmin,wmin],[0 0 0 0],Sigma_mvn))/...
        (mvnpdf([zmax,wmax],[0 0],Sigw)-mvnpdf([zmin,wmin],[0 0],Sigw));
           
end

end

And in .rs file, I define the transition probabilities as, e.g.
! s1_tp_1_2 =1-tp_func2(alpha1,alpha2,rho_v,pev11_0,pev12_1,pev22_0,-inf,b1_0,-inf,b2_0,-inf,tau1,-inf,tau2)

Thank you for your help,
Khoi

jmaih · 2024-10-11T17:59:43Z

Hi Khoi, I guess you are able to evaluate your function directly in Matlab. If so, one thing you may to look at is the implication of the steady state of the model for the evaluation of your function. You can always place a keyboard inside your function so that when RISE calls it, you will be able to easily see what inputs the function receives and whether it returns a finite value. Please try this and let me know, J-

…

On Fri, Oct 11, 2024 at 4:52 PM khoilx ***@***.***> wrote: Hello Prof. Maih, Yes, I pass the alien function correctly (I declared as: m=rise('cmr_mod','steady_state_file','model_ss',... 'rise_flags',model_settings,... 'alien_list','tp_func2');). I'm sending you the 'tp_func2.m' function file. Could you please check it for me? function y=tp_func2(z1,z2,z3,a1,a2,a3,xmin,xmax,ymin,ymax,zmin,zmax,wmin,wmax,varargin) if ~isreal(z1) && ~isreal(z2) && ~isreal(z3) iz1=imag(z1); iz2 = imag(z2); iz3 = imag(z3); z1=real(z1); z2=real(z2); z3=real(z3); if abs(iz1)>1e-12 && abs(iz2)>1e-12 && abs(iz3)>1e-12 y=nan; return end end Pev_sq = [a1 a2; a2 a3]; Pvv = [1, z3; z3, 1]; A = [z1 0; 0 z2]; try % Try solving Sigw = A*Sigw*A'+Pvv for Sigw. Sigw = dlyap(A,Pvv); catch % If A is not stable, throw error and skip all calculations. warning('A is not stable. Assigning 0 to likelihood'); Sigw = eye(2)*.1e10; %replace s.Sigw by garbage %lik = 0; %replace s.lik by 0 end %compute the covariance Pvve Pvve = Pvv - Pev_sq; %compute covariance matrix for 4d multivariate normal in the calculation of time-varying transition Sigma_mvn = [Sigw, Sigw*A'; A*Sigw, Pvve + A*Sigw*A']; if isempty(varargin) y=main(); else if length(varargin)==2 if ~strcmp(varargin{1},'diff') error('Wrong call to function') end switch varargin{2} case 1 y=linhbin_z(); case {2,3,4,5} error('Higher-order derivatives not implemented') end else error('Wrong call to function') end end function y=main() y=mvncdf([xmin,ymin,zmin,wmin],[xmax,ymax,zmax,wmax],[0 0 0 0],Sigma_mvn)/mvncdf([zmin,wmin],[zmax,wmax],[0 0],Sigw); end function fz=linhbin_z() fz=(mvnpdf([xmax,ymax,zmax,wmax],[0 0 0 0],Sigma_mvn)-mvnpdf([xmin,ymin,zmin,wmin],[0 0 0 0],Sigma_mvn))/... (mvnpdf([zmax,wmax],[0 0],Sigw)-mvnpdf([zmin,wmin],[0 0],Sigw)); end end And in .rs file, I define the transition probabilities as, e.g. ! s1_tp_1_2 =1-tp_func2(alpha1,alpha2,rho_v,pev11_0,pev12_1,pev22_0,-inf,b1_0,-inf,b2_0,-inf,tau1,-inf,tau2) Thank you for your help, Khoi — Reply to this email directly, view it on GitHub <#189 (comment)>, or unsubscribe <https://github.com/notifications/unsubscribe-auth/AATKBT6BGFXVEZOB4GKN7RDZ27Q2LAVCNFSM6AAAAABPX2BP4GVHI2DSMVQWIX3LMV43OSLTON2WKQ3PNVWWK3TUHMZDIMBXGU4DCOBWGQ> . You are receiving this because you commented.Message ID: ***@***.***>

khoilx · 2024-10-12T03:21:29Z

Hello Prof. Maih,
I realize that I can only solve the model if I use the .m function with only scalar (the normpdf of a value) to calculate the transitional probabilities, it doesn't work (as retcode shows "3" (Problem in Transition Matrix)) when I declare command of multiple variables cdf (mvdcdf). Do you have any suggestion?

jmaih · 2024-10-12T13:35:30Z

Hi Khoi, Whatever you do inside the alien function is irrelevant to RISE as long as it returns a scalar and the scalar is finite. Is this the case?

…

On Sat, Oct 12, 2024 at 5:21 AM khoilx ***@***.***> wrote: Hello Prof. Maih, I realize that I can only solve the model if I use the .m function with only scalar (the normpdf of a value) to calculate the transitional probabilities, it doesn't work (as recode shows "3" (Problem in Transition Matrix)) when I declare command of multiple variables cdf (mvdcdf). Do you have any suggestion? — Reply to this email directly, view it on GitHub <#189 (comment)>, or unsubscribe <https://github.com/notifications/unsubscribe-auth/AATKBT2WBDEJM22HULFLA4TZ3CIU7AVCNFSM6AAAAABPX2BP4GVHI2DSMVQWIX3LMV43OSLTON2WKQ3PNVWWK3TUHMZDIMBYGMZDSNBZGU> . You are receiving this because you commented.Message ID: ***@***.***>

khoilx · 2024-10-12T14:51:56Z

So what is the way to use the "mvncdf" to calculate the cdf of a vector in the alien function? The output of "mvncdf" is scalar, but its input needs mean vector and covariance matrix.

jmaih · 2024-10-12T15:03:01Z

RISE is not going to give you scalar inputs and expects that you return a scalar output. You are free to list the covariance elements as inputs but RISE doesn't need to know how those inputs are going to be used.

…

On Sat, Oct 12, 2024, 10:52 khoilx ***@***.***> wrote: So what is the way to use the "mvncdf" to calculate the cdf of a vector in the alien function? The output of "mvncdf" is scalar, but its input needs mean vector and covariance matrix. — Reply to this email directly, view it on GitHub <#189 (comment)>, or unsubscribe <https://github.com/notifications/unsubscribe-auth/AATKBT6LXTE6MCAUHDFRHRDZ3EZSFAVCNFSM6AAAAABPX2BP4GVHI2DSMVQWIX3LMV43OSLTON2WKQ3PNVWWK3TUHMZDIMBYGU4TAOJVG4> . You are receiving this because you commented.Message ID: ***@***.***>

jmaih · 2024-10-12T15:03:57Z

RISE will give you scalar inputs, I meant to write.

…

On Sat, Oct 12, 2024, 11:02 Junior Maih ***@***.***> wrote: RISE is not going to give you scalar inputs and expects that you return a scalar output. You are free to list the covariance elements as inputs but RISE doesn't need to know how those inputs are going to be used. On Sat, Oct 12, 2024, 10:52 khoilx ***@***.***> wrote: > So what is the way to use the "mvncdf" to calculate the cdf of a vector > in the alien function? The output of "mvncdf" is scalar, but its input > needs mean vector and covariance matrix. > > — > Reply to this email directly, view it on GitHub > <#189 (comment)>, > or unsubscribe > <https://github.com/notifications/unsubscribe-auth/AATKBT6LXTE6MCAUHDFRHRDZ3EZSFAVCNFSM6AAAAABPX2BP4GVHI2DSMVQWIX3LMV43OSLTON2WKQ3PNVWWK3TUHMZDIMBYGU4TAOJVG4> > . > You are receiving this because you commented.Message ID: > ***@***.***> >

khoilx · 2024-10-12T15:19:46Z

Dear Prof. Maih,
Yes, I used scalar inputs and arrange it into matrix in the alien function to calculate the multivariate normal distribution value.
Could you please check if my alien function is declared correctly?

`function y=tp_func2(z1,z2,z3,a1,a2,a3,xmin,xmax,ymin,ymax,zmin,zmax,wmin,wmax,varargin)

if ~isreal(z1) && ~isreal(z2) && ~isreal(z3)

iz1=imag(z1);

iz2 = imag(z2);

iz3 = imag(z3);

z1=real(z1);

z2=real(z2);

z3=real(z3);

if abs(iz1)>1e-12 && abs(iz2)>1e-12 && abs(iz3)>1e-12

y=nan;

return

end
end

Pev_sq = [a1 a2; a2 a3];

Pvv = [1, z3; z3, 1];

A = [z1 0; 0 z2];

try
% Try solving Sigw = ASigwA'+Pvv for Sigw.
Sigw = dlyap(A,Pvv);
catch
% If A is not stable, throw error and skip all calculations.
warning('A is not stable. Assigning 0 to likelihood');
Sigw = eye(2)*.1e10; %replace s.Sigw by garbage
%lik = 0; %replace s.lik by 0
end

%compute the covariance Pvve
Pvve = Pvv - Pev_sq;

%compute covariance matrix for 4d multivariate normal in the calculation of time-varying transition
Sigma_mvn = [Sigw, SigwA'; ASigw, Pvve + ASigwA'];

if isempty(varargin)

y=main();
else

if length(varargin)==2

if ~strcmp(varargin{1},'diff')
    
    error('Wrong call to function')
    
end

switch varargin{2}
    
    case 1
        
        y=linhbin_z();
        
    case {2,3,4,5}
        
        error('Higher-order derivatives not implemented')
        
end

else

error('Wrong call to function')

end
end

function y=main()

y=mvncdf([xmin,ymin,zmin,wmin],[xmax,ymax,zmax,wmax],[0 0 0 0],Sigma_mvn)/mvncdf([zmin,wmin],[zmax,wmax],[0 0],Sigw);

end

function fz=linhbin_z()

fz=(mvnpdf([xmax,ymax,zmax,wmax],[0 0 0 0],Sigma_mvn)-mvnpdf([xmin,ymin,zmin,wmin],[0 0 0 0],Sigma_mvn))/...
    (mvnpdf([zmax,wmax],[0 0],Sigw)-mvnpdf([zmin,wmin],[0 0],Sigw));

end`
And the transitional probability is defined in .rs file as following:
! s1_tp_1_2 =1-tp_func2(alpha1,alpha2,rho_v,pev11_0,pev12_1,pev22_0,-inf,b1_0,-inf,b2_0,-inf,tau1,-inf,tau2)

khoilx · 2024-10-12T15:26:04Z

The function uses the scalar input (e.g. 1-tp_func2(0.5, 0.5, 0.5, 0.5, 0.5, 0.5, -inf, 0, -inf, 0, -inf, 0, -inf, 0)), and give a scalar result (e.g. 0.25).

jmaih · 2024-10-12T15:52:53Z

1. if 0 and -inf are not going to change, there is no point in writing them as inputs. You can just define them locally inside the function 2. It is not about whether the function runs correctly when you call it with the inputs 0.5, 0.5, .... 3. What matters is when RISE calls the function, particularly the inputs RISE calls the function with. It could well be the case that RISE calls the function with extreme values such that the returned output is not finite. 4. This is why I suggested that you run RISE putting a break inside the function so that you can inspect the inputs received and check whether the output computed is finite and is between 0 and 1

…

On Sat, Oct 12, 2024 at 5:26 PM khoilx ***@***.***> wrote: The function uses the scalar input (e.g. 1-tp_func2(0.5, 0.5, 0.5, 0.5, 0.5, 0.5, -inf, 0, -inf, 0, -inf, 0, -inf, 0)), and give a scalar result (e.g. 0.25). — Reply to this email directly, view it on GitHub <#189 (comment)>, or unsubscribe <https://github.com/notifications/unsubscribe-auth/AATKBT6UQLQPY3TIZRHAZHLZ3E5SHAVCNFSM6AAAAABPX2BP4GVHI2DSMVQWIX3LMV43OSLTON2WKQ3PNVWWK3TUHMZDIMBYGYYDCMBTGE> . You are receiving this because you commented.Message ID: ***@***.***>

Provide feedback

Saved searches

Use saved searches to filter your results more quickly

Problem in transition matrix #189

Problem in transition matrix #189

khoilx commented Oct 10, 2024

jmaih commented Oct 11, 2024 via email

khoilx commented Oct 11, 2024

jmaih commented Oct 11, 2024 via email

khoilx commented Oct 12, 2024 •

edited

Loading

jmaih commented Oct 12, 2024 via email

khoilx commented Oct 12, 2024

jmaih commented Oct 12, 2024 via email

jmaih commented Oct 12, 2024 via email

khoilx commented Oct 12, 2024

khoilx commented Oct 12, 2024

jmaih commented Oct 12, 2024 via email

Problem in transition matrix #189

Problem in transition matrix #189

Comments

khoilx commented Oct 10, 2024

jmaih commented Oct 11, 2024 via email

khoilx commented Oct 11, 2024

jmaih commented Oct 11, 2024 via email

khoilx commented Oct 12, 2024 • edited Loading

jmaih commented Oct 12, 2024 via email

khoilx commented Oct 12, 2024

jmaih commented Oct 12, 2024 via email

jmaih commented Oct 12, 2024 via email

khoilx commented Oct 12, 2024

khoilx commented Oct 12, 2024

jmaih commented Oct 12, 2024 via email

khoilx commented Oct 12, 2024 •

edited

Loading