Home > wafo > wstats > wtinv.m

wtinv

PURPOSE ^

Inverse of the Student's T distribution function

SYNOPSIS ^

x = wtinv(F,df)

DESCRIPTION ^

 WTINV Inverse of the Student's T distribution function 
  
  CALL:  x = wtinv(F,df); 
  
    x   = inverse cdf for Student's T distribution evaluated at F. 
    df  = degrees of freedom 
     
   Example: 
     F = linspace(0,1,100); 
     x = wtinv(F,1); 
     plot(F,x)

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SOURCE CODE ^

001 function x = wtinv(F,df) 
002 %WTINV Inverse of the Student's T distribution function 
003 % 
004 % CALL:  x = wtinv(F,df); 
005 % 
006 %   x   = inverse cdf for Student's T distribution evaluated at F. 
007 %   df  = degrees of freedom 
008 %    
009 %  Example: 
010 %    F = linspace(0,1,100); 
011 %    x = wtinv(F,1); 
012 %    plot(F,x) 
013  
014  
015 % tested on matlab 5.3 
016 %History: 
017 %revised pab 22.05.2003 
018 % - added new method for df==2 
019 %revised pab 29.10.2000 
020 % adapted from stixbox 
021 % -added nargchk, comnsize 
022 % - added check on F, df + some code from wnormplot when df>mxdf 
023 %       Anders Holtsberg, 18-11-93 
024 %       Copyright (c) Anders Holtsberg 
025  
026 error(nargchk(2,2,nargin)) 
027 [errorcode F,df] = comnsize(F,df); 
028 if errorcode>0, 
029   error('F and df must be of common size or scalar'); 
030 end 
031 x = zeros(size(F)); 
032  
033 mxdf = 10^2; 
034  
035 ok = (df>0 & floor(df)==df & 0 <=F & F<=1  ); 
036 region1 = ((0 < F) & (F < 1)); 
037  
038 k1 = find((df == 1) & region1 & ok); 
039 if any(k1) 
040   x(k1) = tan(pi * (F(k1) - 0.5)); 
041 end 
042  
043 k2 = find((df == 2)& region1 & ok); 
044 if (any(k2)) 
045   R = (2*F(k2)-1); 
046   x(k2) = sqrt(2)*R./sqrt(1-R.*R); 
047 end 
048  
049 k = find(ok & region1 &  2 < df  & df<=mxdf); 
050 if any(k), 
051   s = (F(k)<0.5);  
052   tmp = F(k) + (1-2*F(k)).*s; 
053   tmp2 = wbetainv(1-(2*(1-tmp)),1/2,df(k)/2); 
054   tmp3 = tmp2.*df(k)./((1-tmp2)); 
055   x(k) = (1-2*s).*sqrt(tmp3); 
056 end 
057 k=find(ok & region1 & df>mxdf); 
058 if any(k), % pab 01.11.2000, added from wnormplot  
059   x(k) = wnorminv(F(k)); 
060   k0=find(df(k)<inf); 
061   if any(k0) 
062     k1=k(k0); 
063     Y=x(k1); 
064     g1=1/4*(Y.^3+Y); 
065     g2=1/96*(5*Y.^5+16*Y.^3+3*Y); 
066     g3=1/384*(3*Y.^7+19*Y.^5+17*Y.^3-15*Y); 
067     g4=1/92160*(79*Y.^9+776*Y.^7+1482*Y.^5-1920*Y.^3-945*Y); 
068     x(k1)=Y+g1./df(k1)+g2./df(k1).^2+g3./df(k1).^3+g4./df(k1).^4; 
069   end 
070 end 
071   
072  
073 k2=find(ok&F==0); 
074 if any(k2) 
075   x(k2)=-inf; 
076 end 
077  
078 k3=find(ok&F==1); 
079 if any(k3) 
080   x(k3)=inf; 
081 end 
082  
083  
084  
085  
086 k4=find(~ok); 
087 if any(k2) 
088   x(k4)=NaN; 
089 end 
090  
091

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