Home > wafo > wstats > weib2dcdf.m

weib2dcdf

PURPOSE ^

Joint 2D Weibull cumulative distribution function

SYNOPSIS ^

[y ,tol1] = weib2dcdf(x1,x2,param,condon,tol)

DESCRIPTION ^

 WEIB2DCDF Joint 2D Weibull cumulative distribution function   
  
   CALL: [F tol]= weib2dcdf(x1,x2,param,condon,rtol) 
  
      F    = joint CDF, i.e., Prob(X1<x1 X2<x2) 
     tol   = estimated absolute tolerance 
     x1,x2 = coordinates 
     param = [ A1, B1,A2,B2,C12] the parameters of the distribution 
    condon = 0 it returns the the regular cdf of X1 and X2 (default) 
             1 it returns the conditional cdf given X1 
             2 it returns the conditional cdf given X2 
      rtol = specified relative tolerance (Default 1e-3)  
  
  The size of F is the common size of X1 and X2. 
  The CDF is defined by its PDF: 
  
   f(X1,X2)=B1*B2*xn1^(B1-1)*xn2^(B2-1)/A1/B1/N*... 
            exp{-[xn1^B1 +xn2^B2 ]/N }*I0(2*C12*xn1^(B1/2)*xn2^(B2/2)/N)  
    where  
       N=1-C12^2, xn1=X1/A1,  xn2=X2/A2 and  
       I0 is the modified  bessel function of zeroth order. 
  
   (The marginal distribution of X1 is weibull with parameters A1 and B1)  
   (The marginal distribution of X2 is weibull with parameters A2 and B2)  
  
  Example: 
    x = linspace(0,6,200); [X1,X2] = meshgrid(x);  
    phat = [2 2  3 2.5 .8]; 
    f = weib2dpdf(X1,X2,phat); 
    contour(x,x,f), hold on, 
    plot( [0  2  2 0 0], [0 0 1 1 0],'g-'), hold off % Mark the region 
    weib2dcdf(2,1,phat)  % Calculate the probability of marked region 
   
  See also  weib2dpdf, wweibpdf, quad2dg, gaussq

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SOURCE CODE ^

001 function [y ,tol1] = weib2dcdf(x1,x2,param,condon,tol) 
002 %WEIB2DCDF Joint 2D Weibull cumulative distribution function   
003 % 
004 %  CALL: [F tol]= weib2dcdf(x1,x2,param,condon,rtol) 
005 % 
006 %     F    = joint CDF, i.e., Prob(X1<x1 X2<x2) 
007 %    tol   = estimated absolute tolerance 
008 %    x1,x2 = coordinates 
009 %    param = [ A1, B1,A2,B2,C12] the parameters of the distribution 
010 %   condon = 0 it returns the the regular cdf of X1 and X2 (default) 
011 %            1 it returns the conditional cdf given X1 
012 %            2 it returns the conditional cdf given X2 
013 %     rtol = specified relative tolerance (Default 1e-3)  
014 % 
015 % The size of F is the common size of X1 and X2. 
016 % The CDF is defined by its PDF: 
017 % 
018 %  f(X1,X2)=B1*B2*xn1^(B1-1)*xn2^(B2-1)/A1/B1/N*... 
019 %           exp{-[xn1^B1 +xn2^B2 ]/N }*I0(2*C12*xn1^(B1/2)*xn2^(B2/2)/N)  
020 %   where  
021 %      N=1-C12^2, xn1=X1/A1,  xn2=X2/A2 and  
022 %      I0 is the modified  bessel function of zeroth order. 
023 % 
024 %  (The marginal distribution of X1 is weibull with parameters A1 and B1)  
025 %  (The marginal distribution of X2 is weibull with parameters A2 and B2)  
026 % 
027 % Example: 
028 %   x = linspace(0,6,200); [X1,X2] = meshgrid(x);  
029 %   phat = [2 2  3 2.5 .8]; 
030 %   f = weib2dpdf(X1,X2,phat); 
031 %   contour(x,x,f), hold on, 
032 %   plot( [0  2  2 0 0], [0 0 1 1 0],'g-'), hold off % Mark the region 
033 %   weib2dcdf(2,1,phat)  % Calculate the probability of marked region 
034 %  
035 % See also  weib2dpdf, wweibpdf, quad2dg, gaussq 
036  
037 % tested on: matlab5.1 
038 %history: 
039 % revised pab Dec2003 
040 %  call weibstat -> wweibstat 
041 % revised pab 29.10.2000 
042 % - updated to wstats 
043 % - added example 
044 %  Per A. Brodtkorb 17.11.98 
045  
046 % NB! weibpdf must be modified to correspond to 
047 % pdf=x^(b-1)/a^b*exp(-(x/a)^b)  
048  
049  
050 error(nargchk(3,5,nargin)) 
051 if (nargin <5)| isempty(tol),   tol=1e-3; end 
052 if (nargin <4)| isempty(condon),  condon=0;end 
053 if nargin < 3,  
054   error('Requires 3 input arguments.');  
055 elseif prod(size((param)))~=5|isempty(param), 
056    error('input argument 3 must have 5 entries.');  
057 else 
058   a1=param(1);  b1=param(2);  
059   a2=param(3);  b2=param(4);    c12=param(5);  
060 end 
061  
062 [errorcode x1 x2 ] = comnsize(x1,x2); 
063 if  errorcode > 0 
064   error('Requires non-scalar arguments to match in size.'); 
065 end 
066 y = zeros(size(x1)); 
067 tol1=y; 
068  
069 if 0 
070   % This is a trick to get the html documentation correct. 
071   k = weib2dpdf(1,1,2,3); 
072 end 
073  
074 switch condon 
075   case 0,% regular cdf    
076     k = find( x1 > 0 & x2>0); 
077     if 1,  
078        if any(k),  
079      [y(k) tol1(k)]= quad2dg('weib2dpdf',0,x1(k),0 , x2(k),  tol,a1,b1,a2,b2,c12,condon);   
080        end 
081         
082         
083      else% divide the integral in to several parts this is not correct yet 
084       [m1 v1]=wweibstat(a1,b1); [m2 v2]=wweibstat(a2,b2); 
085        x1s=min(m1,x1); x2s=min(m2,x2); 
086        if any(k) 
087      [y(k) tol1(k)]=quad2dg('weib2dpdf',0,x1s(k),0 , x2s(k),tol/2,param,condon); 
088        end 
089        k1=find( x1(k)>m1& x2(k)<=m2); 
090        if any(k1), 
091      [tmp1 tmp2]=quad2dg('weib2dpdf',x1s(k(k1)),  x1(k(k1)),0,   x2s(k(k1)), tol/4,param,condon); 
092      y(k(k1)) =y(k(k1))+tmp1;tol1(k(k1))=tol1(k(k1))+tmp2; 
093        end 
094        k1=find( x1(k)<=m1& x2(k)>m2); 
095        if any(k1), 
096      [tmp1 tmp2]=quad2dg('weib2dpdf',0,x1s(k(k1)), x2s(k(k1)),   x2(k(k1)),   tol/4,param,condon); 
097      y(k(k1)) =y(k(k1)) +tmp1;tol1(k(k1))=tol1(k(k1))+tmp2; 
098        end 
099        k1=find(m1<x1(k)& m2<x2(k)); 
100        if any(k1), 
101      [tmp1 tmp2]=quad2dg('weib2dpdf',x1s(k(k1)),x1(k(k1)), x2s(k(k1)),   x2(k(k1)),   tol/4,param,condon); 
102      y(k(k)) =y(k(k)) +tmp1;tol1(k(k))=tol1(k(k))+tmp2; 
103        end 
104      end 
105    case 1,%conditional CDF given x1   
106      k = find( (x1 > 0) & (x2>0) & (c12(ones(size(x2))) ~=0 )); 
107      if any(k), 
108        [y(k) tol1(k)]=gaussq('weib2dpdf',0,x2(k),tol,[],x1(k),a2,b2,a1,b1,c12,2);%param([3:4 1:2 5]),2);    
109      end 
110      k = find( (x1==0)| (c12(ones(size(x2))) ==0)); 
111      if any(k), 
112        y(k) =wweibcdf(x2(k),param(3),param(4));  
113      end 
114      %for ix=1:length(x1(:)),       
115      %[y(ix) tol1(ix)]=quadg('weib2dpdf',0,x2(ix),tol,[],x1(ix),param([3:4 1:2 5]),2);       
116      %end 
117    case 2,%conditional CDF given x2   
118      k = find( (x1 > 0) & (x2>0) & (c12(ones(size(x2))) ~=0)); 
119      if any(k), 
120        [y(k) tol1(k)]=gaussq('weib2dpdf',0,x1(k),tol,[],x2(k),a1,b1,a2,b2,c12,2);%param,2);  
121      end 
122      k = find( (x2==0)| (c12(ones(size(x2))) ==0)); 
123       if any(k), 
124        y(k) =wweibcdf(x1(k),param(1),param(2));  
125      end 
126      %for ix=1:length(x2(:)),  
127       % [y(ix) tol1(ix)]=quadg('weib2dpdf',0,x1(ix),tol,[],x2(ix),param,2);       
128      %end 
129  end     
130   
131 % make sure 0 <= y<=1  
132 k2=find(y<0); 
133 if any(k2) 
134   y(k2)=zeros(size(k2)); 
135 end 
136 k2=find(y>1); 
137 if any(k2) 
138   y(k2)=ones(size(k2)); 
139 end 
140 if any(isnan(y)), 
141   disp(['Warning: there are : ', num2str(sum(isnan(y))),' NaNs']); 
142 end 
143

Mathematical Statistics
Centre for Mathematical Sciences
Lund University with Lund Institute of Technology

Comments or corrections to the WAFO group


Generated on Thu 06-Oct-2005 02:21:16 for WAFO by m2html © 2003