MDIST2DPDF Joint 2D PDF due to Plackett given as f{x1}*f{x2}*G(x1,x2;Psi). CALL: f = mdist2dpdf(x1,x2,phat,condon) f = PDF evalutated at points (x1 , x2) phat = parameter structure containing x{1} = Phat1 marginal parameters of X1 x{2} = Phat2 marginal parameters of X2 x{3} = Psi the interaction parameter between x1 and x2. dist = list of marginal distributions of x1 and x2, respectively Options are: 'tgumbel', 'gumbel', 'lognormal','rayleigh','weibull','gamma'. condon = 0 regular pdf is returned (default) 1 conditional pdf of H given V is returned 2 conditional pdf of V given H is returned Example: 2D Weibull Rayleigh with parameters [2 3] and 3, respectively for the marginal distributions, and a interaction parameter of 10: phat.x={[2 3],3,10}; phat.dist={'weibull','rayleigh'}; f = mdist2dpdf(3,4,phat); See also mdist2dcdf, mdist2dfit, mdist2drnd
Check if all input arguments are either scalar or of common size. | |
Gamma cumulative distribution function | |
Gamma probability density function | |
Gumbel cumulative distribution function. | |
Gumbel probability density function. | |
Lognormal cumulative distribution function | |
Lognormal probability density function | |
Rayleigh cumulative distribution function | |
Rayleigh probability density function | |
Weibull cumulative distribution function | |
Weibull probability density function | |
Display message and abort function. | |
Convert string to lowercase. | |
Compare strings. |
Inverse of the conditional cdf of X2 given X1. | |
MDIST log-likelihood function. | |
Joint 2D PDF due to Plackett given as f{x1}*f{x2}*G(x1,x2;Psi). | |
Mean and variance for the MDIST2D distribution. |
001 function y = mdist2dpdf(V,H,phat,condon) 002 %MDIST2DPDF Joint 2D PDF due to Plackett given as f{x1}*f{x2}*G(x1,x2;Psi). 003 % 004 % CALL: f = mdist2dpdf(x1,x2,phat,condon) 005 % 006 % f = PDF evalutated at points (x1 , x2) 007 % phat = parameter structure containing 008 % x{1} = Phat1 marginal parameters of X1 009 % x{2} = Phat2 marginal parameters of X2 010 % x{3} = Psi the interaction parameter between x1 and x2. 011 % dist = list of marginal distributions of x1 and x2, respectively 012 % Options are: 'tgumbel', 'gumbel', 013 % 'lognormal','rayleigh','weibull','gamma'. 014 % condon = 0 regular pdf is returned (default) 015 % 1 conditional pdf of H given V is returned 016 % 2 conditional pdf of V given H is returned 017 % 018 % Example: 2D Weibull Rayleigh with parameters [2 3] and 3, 019 % respectively for the marginal distributions, and a interaction 020 % parameter of 10: 021 % 022 % phat.x={[2 3],3,10}; 023 % phat.dist={'weibull','rayleigh'}; 024 % f = mdist2dpdf(3,4,phat); 025 % 026 % See also mdist2dcdf, mdist2dfit, mdist2drnd 027 028 029 % References: 030 % Plackett, R. L. (1965) "A class of bivariate distributions." 031 % J. Am. Stat. Assoc. 60. 516-22 032 % [1] Michel K. Ochi, 033 % OCEAN TECHNOLOGY series 6 034 % "OCEAN WAVES, The stochastic approach", Cambridge 035 % 1998 pp. 133-134. 036 037 038 % tested on: matlab 5.2 039 % history 040 % revised pab 20.10.2000 041 % - updated to new wstats toolbox 042 % revised pab 8.11.1999 043 % - updated header info 044 % - changed phat from vectro to structure 045 % Per A. Brodtkorb 28.01.99 046 047 048 error(nargchk(3,4,nargin)) 049 if nargin <4 |isempty(condon), condon =0;end 050 [errorcode V H ] = comnsize(V,H); 051 if errorcode > 0 052 error ('x1 and x2 must be of common size or scalar'); 053 end 054 055 VDIST=lower(phat.dist{1}); 056 HDIST=lower(phat.dist{2}); 057 058 059 psi=phat.x{3}; % interaction parameter 060 PV=phat.x{1}; 061 PH=phat.x{2}; 062 063 064 065 y = zeros(size(V)); 066 067 068 if strcmp('gu', VDIST(1:2)), 069 if strcmp('gu',HDIST(1:2)), 070 k=find(H>-inf&V>-inf); 071 else 072 k=find(H>0&V>-inf); 073 end 074 elseif strcmp('gu',HDIST(1:2)), 075 k = find(H>-inf & V > 0 ); 076 else 077 k = find( H>0 &V > 0); 078 end 079 080 if any(k), 081 Fv=dist1dcdffun(V(k),PV, VDIST(1:2)); 082 fv=dist1dpdffun(V(k),PV, VDIST(1:2)); 083 Fh=dist1dcdffun(H(k),PH, HDIST(1:2)); 084 fh=dist1dpdffun(H(k),PH, HDIST(1:2)); 085 086 tmp=1+(Fv+Fh).*(psi-1); 087 y(k)=psi.*((psi-1).*(Fv+Fh-2.*Fv.*Fh)+1)./(sqrt(tmp.^2-4.*psi.*(psi-1).*Fv.*Fh).^3); 088 switch condon 089 case 0, y(k)=y(k).*fv.*fh; 090 case 1, y(k)=y(k).*fh; 091 case 2, y(k)=y(k).*fv; 092 case 3, % secret option used by mdist2dstat: returns v*f(v|h) 093 y(k)=V(k).*y(k).*fv; 094 case 4, % secret option used by weib2dstat: returns v^2*f(v|h) 095 y(k)=V(k).^2.*y(k).*fv; 096 end 097 %plot(V(k),y(k),'.') 098 end 099 100 function cdf1=dist1dcdffun(H,Ah,dist2 ) 101 switch dist2(1:2) 102 case 'ra', cdf1=wraylcdf(H,Ah); 103 case 'we' , cdf1=wweibcdf(H,Ah(1),Ah(2)); 104 case 'gu' , cdf1=wgumbcdf(H,Ah(1),Ah(2),0); 105 case 'tg' , cdf1=wgumbcdf(H,Ah(1),Ah(2),1); 106 case 'ga' , cdf1=wgamcdf(H,Ah(1),Ah(2)); 107 case 'lo' , cdf1=wlogncdf(H,Ah(1),Ah(2)); 108 otherwise, error('unknown distribution') 109 end 110 return 111 function pdf1=dist1dpdffun(H,Ah,dist2 ) 112 switch dist2(1:2) 113 case 'ra', pdf1=wraylpdf(H,Ah); 114 case 'we' , pdf1=wweibpdf(H,Ah(1),Ah(2)); 115 case 'gu' , pdf1=wgumbpdf(H,Ah(1),Ah(2),0); 116 case 'tg' , pdf1=wgumbpdf(H,Ah(1),Ah(2),1); 117 case 'ga' , pdf1=wgampdf(H,Ah(1),Ah(2)); 118 case 'lo' , pdf1=wlognpdf(H,Ah(1),Ah(2)); 119 otherwise, error('unknown distribution') 120 end 121 return 122 123 124 125 126 127 128
Comments or corrections to the WAFO group