Contents
function stat = stats_estimator_A(n,k,varargin)
% stat = stats_estimator_A(n,k,p) % Gives some statistical properties of the estimator_A % % Inputs % * n : total number of particles % * k : number of survivors % % Parameters % * alpha : confidence interval at alpha % % * P : probability to be estimated % OR % * iter : number of iterations % * k_prime : number of particles above the last threshold % % Outputs % stat : a cell containing % * rel_bias : relative bias % * rel_var : relative variance % * interval : Interval of confidence % * NbIter : number of iterations % % Comments % - % % Authors % * Implementation: Teddy Furon % * Science: Frederic Cerou, Pierre Del Moral, Arnaud Guyader, Teddy Furon % % Contact: teddy.furon@inria.fr % % Copyright INRIA 2010 % % Licence % This software and its subprograms are governed by the CeCILL license % under French law and abiding by the rules of distribution of free software. % See http://www.cecill.info/licences.en.html
Parsing the arguments
error(nargoutchk(0,3,nargout)); % Create an instance of the inputParser class. p = inputParser; % Define inputs that one must pass on every call: p.addRequired('n', @(x)validateattributes(x,{'numeric'},{'integer', 'positive'})); p.addRequired('k', @(x)validateattributes(x,{'numeric'},{'integer', 'positive'})); p.addParamValue('k_prime',1,@(x)validateattributes(x,{'numeric'},{'integer', 'positive'})); p.addParamValue('iter',100,@(x)validateattributes(x,{'numeric'},{'integer', 'positive'})); p.addParamValue('P',10^-7,@(x)validateattributes(x,{'numeric'},{'positive','scalar'})); p.addParamValue('alpha',0.95,@(x)validateattributes(x,{'numeric'},{'positive','scalar'})); % Parse and validate all input arguments. p.StructExpand = true; p.KeepUnmatched = true; p.parse(n,k,varargin{:}); defaulted = p.UsingDefaults; h = p.Results;
Input argument "n" is undefined. Error in ==> stats_estimator_A at 57 p.parse(n,k,varargin{:});
Which mode?
if nnz(strcmp('P',defaulted))>0
Operational mode
'P' was not passed a real estimator was runned before
p0 = h.k/h.n; n0 = h.iter-1; r0 = h.k_prime/h.n; Prob = p0^n0*r0;
else
Theoretical mode
'P' was passed We want an estimate of the properties for such a probability to be estimated
p0 = h.k/h.n; n0 = floor(log(h.P)/log(p0)); r0 = h.P*p0^(-n0); Prob = h.P;
end
Statistical properties
stat.NbIter = n0+1; % relative bias of the estimation stat.rel_bias = n0*(1-p0)/p0/h.n; % relative variance of the estimation stat.rel_var = (n0*(1-p0)/p0 + (1-r0)/r0)/h.n; % Interval of alpha% confidence La = sqrt(2)*erfinv(h.alpha); stat.interval = Prob*(1-stat.rel_bias+sqrt(stat.rel_var)*[-La,La]); stat.NbIter = n0+1;