function stat_ds=cosmo_phase_stat(ds,varargin)
% Compute phase perturbation, or opposition sum or product phase statistic
%
% stat_ds=cosmo_phase_stat(ds,...)
%
% Inputs:
% ds dataset struct with fields:
% .samples PxQ complex matrix for P samples (trials,
% observations) and Q features (e.g. combinations
% of time points, frequencies and channels)
% .sa.targets Px1 array with trial conditions.
% There must be exactly two conditions, thus
% .sa.targets must have exactly two unique
% values. A balanced number of samples is
% requires, i.e. each of the two unique values in
% .sa.targets must occur equally often.
% .sa.chunks Px1 array indicating which samples can be
% considered to be independent. It is required
% that all samples are independent, therefore
% all values in .sa.chunks must be different from
% each other
% .fa } optional feature attributes
% .a } optional sample attributes
% 'output',p Return statistic, one of the following:
% - 'pbi': phase bifurcation index
% - 'pos': phase opposition sum
% - 'pop': phase opposition product
% 'samples_are_unit_length',u (optional)
% If u==true, then all elements in ds.samples
% are assumed to be already of unit length. If
% this is indeed true, this can speed up the
% computation of the output.
% 'check_dataset',c (optional, default=true)
% if c==false, there is no check for consistency
% of the ds input.
%
% Output:
% stat_ds struct with fields
% .samples 1xQ array with 'pbi', 'pos', or 'pop' function
% .a } if present in the input, then the output
% .fa } contains these fields as well
%
% Notes:
% - if a dataset is not balanced for number of trials, consider using
% cosmo_balance_dataset to balance it.
%
% See also: cosmo_balance_dataset, cosmo_phase_itc
%
% # For CoSMoMVPA's copyright information and license terms, #
% # see the COPYING file distributed with CoSMoMVPA. #
