AtkinsonIneq(w, epsilon)

function y = AtkinsonIneq(w, epsilon)

% The function computes the Atkinson Inequality Index for the wealth vector w associated to single individuals.
% w must be non-negative (with at least one i such that w(i) > 0).
% epsilon must be a real in the interval (0, 1)
% The Atkinson Inequality Index is a measure of inequality, i.e. a measure of wealth concentration.
%
% http://en.wikipedia.org/wiki/Lorenz_curve
% http://en.wikipedia.org/wiki/Gini_coefficient
% http://en.wikipedia.org/wiki/Theil_index
% http://en.wikipedia.org/wiki/Atkinson_index
%
% -*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*
% -*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*
%
% % Example 1: Uniform U(0, 1)
% N = 1000;                                           % Number of individuals
% w = rand(N, 1);                                     % Wealth extracted from a Uniform U(0, 1)
% epsilon = rand;
% y = AtkinsonIneq(w, epsilon)
% y =
% 
%     0.1183
%
% -*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*
% -*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*
%
% % Example 2: Standard Normal z(0, 1)
% N = 1000;                                           % Number of individuals
% mu = 4;
% w = mu + randn(N, 1);                               % Wealth extracted from a Normal N(mu, 1)
% w = abs(w);                                         % Be careful: all values must be strictly positive!!!
% epsilon = rand;
% y = AtkinsonIneq(w, epsilon)
% y =
% 
%     0.0198
%
% -*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*
% -*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*
%
% % Example 3: LogNormal logN(0, 1)
% N = 1000;                                           % Number of individuals
% w = exp(randn(N, 1));                               % Wealth extracted from a LogNormal logN(0, 1)
% epsilon = rand;
% y = AtkinsonIneq(w, epsilon)
% y =
% 
%     0.2111
%
% -*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*
% -*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*
%
% % Example 4: Power Law PL(kappaPL, alphaPL)
% N = 1000;                                           % Number of individuals
% kappaPL = 1;
% alphaPL = 2;
% w = kappaPL * (rand(N, 1) .^ (-1/alphaPL));         % Wealth extracted from a Power Law PL(kappaPL, alphaPL)
% epsilon = rand;
% y = AtkinsonIneq(w, epsilon)
% y =
% 
%     0.1022
%
% -*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*
% -*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*
%
% % Example 5: Mixture of LogNormal logN(0, 1) and Power Law PL(kappaPL, alphaPL)
% N = 1000;                                           % Number of individuals
% alpha = floor(0.92 * N);                            % Number of individuals whose wealth is extracted from a LogNormal LogN(0, 1)
% w = exp(randn(alpha, 1));                           % Wealth extracted from a LogNormal LogN(0, 1)
% kappaPL = 1;
% alphaPL = 2;
% w((alpha + 1):N) = kappaPL * (rand(N - alpha, 1) .^ (-1/alphaPL));   % Wealth extracted from a Power Law PL(kappaPL, alphaPL)
% epsilon = rand;
% y = AtkinsonIneq(w, epsilon)
% y =
% 
%     0.2042
%
% -*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*
% -*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*
%
% Written by
% Francesco Pozzi
% 19 May 2008
%

ctrl = isvector(w) & isnumeric(w) & isreal(w) & isscalar(epsilon) & isnumeric(epsilon) & isreal(epsilon);
if ctrl
  w = w(:);
  w = w(~isnan(w) & ~isinf(w));
else
  error('Epsilon and/or Wealth Values incorrect!')
end

ctrl = all(w >= 0) & any(w > 0) & (size(epsilon) == [1, 1]) & (epsilon > 0) & (epsilon <= 1);
if ~ctrl
  error('Epsilon and/or Wealth Values incorrect!')
end

mu = mean(w);
N = length(w);

if epsilon == 1;
  y = 1 - (prod(w) ^ (1 / N)) / mu;
else
  y = 1 - ((sum(w .^ (1 - epsilon)) / N) ^ (1 / (1 - epsilon))) / mu;
end