Home > wafo > wstats > wgpdpdf.m

wgpdpdf

PURPOSE ^

Generalized Pareto probability density function

SYNOPSIS ^

f = wgpdpdf(x,k,s,m);

DESCRIPTION ^

 WGPDPDF Generalized Pareto probability density function 
  
  CALL:  f = wgpdpdf(x,k,s,m); 
   
         f = density function evaluated at x 
         k = shape parameter in the GPD 
         s = scale parameter in the GPD    (default 1) 
         m = location parameter in the GPD (default 0)   
   
  The Generalized Pareto distribution is defined by its cdf 
  
                 1 - (1-k(x-m)/s)^1/k,  k~=0 
   F(x;k,s) = 
                 1 - exp(-(x-m)/s),  k==0 
  
  for x>=m (when k<=0) and 0 <= x-m < s/k (when k>0), s>0. 
  
  Example:  
    x = linspace(0,2,200); 
    p1 = wgpdpdf(x,1.25,1); p2 = wgpdpdf(x,1,1); 
    p3 = wgpdpdf(x,0.75,1); p4 = wgpdpdf(x,0.5,1); 
    plot(x,p1,x,p2,x,p3,x,p4)

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SOURCE CODE ^

001 function f = wgpdpdf(x,k,s,m); 
002 %WGPDPDF Generalized Pareto probability density function 
003 % 
004 % CALL:  f = wgpdpdf(x,k,s,m); 
005 %  
006 %        f = density function evaluated at x 
007 %        k = shape parameter in the GPD 
008 %        s = scale parameter in the GPD    (default 1) 
009 %        m = location parameter in the GPD (default 0)   
010 %  
011 % The Generalized Pareto distribution is defined by its cdf 
012 % 
013 %                1 - (1-k(x-m)/s)^1/k,  k~=0 
014 %  F(x;k,s) = 
015 %                1 - exp(-(x-m)/s),  k==0 
016 % 
017 % for x>=m (when k<=0) and 0 <= x-m < s/k (when k>0), s>0. 
018 % 
019 % Example:  
020 %   x = linspace(0,2,200); 
021 %   p1 = wgpdpdf(x,1.25,1); p2 = wgpdpdf(x,1,1); 
022 %   p3 = wgpdpdf(x,0.75,1); p4 = wgpdpdf(x,0.5,1); 
023 %   plot(x,p1,x,p2,x,p3,x,p4) 
024   
025 % References 
026 %  Johnson N.L., Kotz S. and Balakrishnan, N. (1994) 
027 %  Continuous Univariate Distributions, Volume 1. Wiley.  
028  
029 % Tested on: Matlab 5.3 
030 % History:  
031 % Revised pab oct 2005 
032 % - limits m<x<s/k replaced with 0 < x-m < s/k in help header. 
033 % -fixed a bug for k<0  
034 % revised pab June 2005 
035 % fixed bug for k==0 now fixed 
036 % revised jr 14.08.2001 
037 %  - a bug in the last if-statement condition fixed 
038 %    (thanks to D Eddelbuettel <edd@debian.org>) 
039 % revised pab 24.10.2000 
040 %  - added comnsize, nargchk 
041 %  - added extra check on the scale parameter s 
042 %  - added m + default values on m and s 
043 % added ms 14.06.2000 
044  
045 error(nargchk(2,4,nargin)) 
046  
047 if nargin<4|isempty(m), m=0;end 
048 if nargin<3|isempty(s), s=1;end 
049  
050 [errorcode x k s m] = comnsize(x,k,s,m); 
051 if errorcode > 0 
052     error('x, k s and m must be of common size or scalar.'); 
053 end 
054 epsilon=1e-4; % treshold defining k to zero 
055  
056 xm = x-m; 
057 f = zeros(size(x)); 
058  
059 k0 = find((xm>=0) & (abs(k)<=epsilon) & (s>0)); 
060 if any(k0),  
061    f(k0) = exp(-xm(k0)./s(k0))./s(k0); 
062 end 
063  
064 k1=find((xm>=0) & (k.*xm < s) & (abs(k)>epsilon) & (s>0)); 
065 if any(k1), 
066   f(k1) = (1-k(k1).*xm(k1)./s(k1)).^(1./k(k1)-1)./s(k1); 
067 end 
068  
069    
070 k3 = find(s<=0); 
071 if any(k3), 
072    tmp   = NaN; 
073    f(k3) = tmp(ones(size(k3))); 
074 end 
075  
076  
077  
078  
079  
080  
081  
082  
083  
084  
085  
086  
087

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