Jos

This submission has several problems.
1. The help section is insufficient: add a H1 line, a clear description how to use this function and what to expect, and add an example.

2. The for-loop is not pre-allocated, and therefore quite slow with large n

3. There are several faster alternatives, for instance
n = 4 ;
R = rem(floor([0:((2^n)-1)].'* pow2(0:-1:-n+1)),2)
(adapted from DEC2BIN)

It is nice to see some improvments Franciso! A couple of additional comments: if your sor worried about speed, I suggest you take a look at logical indexing, and replace all the zeros from the struct command by one call to zeros.

Another point: the behavior for 2D matrices should be explained, especially when the groups are specified. I think you did not intend that behavior. For instance, take a look if the following produces what you want them to produce:

A = rand(3,4) ; allstats(A), allstats(A,1:4), allstats(A,1:3)

With respect to the help, you can add examples of the expected output, besides giving only the calling syntax.

Below is the outline of a more optimized code:
function R = allstats(DATA, GROUP)
if isempty(ver('stats'))
    error('Function requires STATS toolbox.') ;
end
if nargin==1,
    GROUP = ones(numel(DATA),1) ;
    Ngroups = 1 ;
    UniqueGroups = 1 ;
else
    if numel(DATA) ~= numel(GROUP)
        error('Number of elements should match') ;
    end
    UniqueGroups = unique(GROUP) ;
    Ngroups = numel(UniqueGroups) ;
end
Z = nan(Ngroups,1) ;
R = struct('min',Z,'max',Z,'mean',Z) ; % etce
for k=1:Ngroups,
    Q = GROUP == UniqueGroups(k) ; % logical indexing
    R(k).min = nanmin(DATA(Q)) ;
    R(k).max = nanmax(DATA(Q)) ;
    R(k).mean = nanmean(DATA(Q)) ;
end

Just what I was looking for

I've tried this small experiment:

N = 100000; % sample size
w = 1:100; % quite steep distribution
w = w / sum(w); % distribution (pdf)
Y_std = sqrt(w .* (1 - w)); % Standard deviations associated to each probability

%%%%%%%%%%%%%%%%%

Y1 = randp(w, N, 1); % Try Jos function
Y1_w = histc(Y1, [1:100]); % absolute frequencies table
Y1_w = Y1_w / N; % transform into relative frequencies
Y1_w_std = sqrt(Y1_w .* (1 - Y1_w)); % Empirical standard deviations associated to each probability

%%%%%%%%%%%%%%%%%

Y2 = gDiscrPdfRnd(w, N, 1); % Try Gianluca Dorini's function
Y2_w = histc(Y2, [1:100]); % absolute frequencies table
Y2_w = Y2_w / N; % transform into relative frequencies
Y2_w_std = sqrt(Y2_w .* (1 - Y2_w)); % Empirical standard deviations associated to each probability

%%%%%%%%%%%%%%%%%

% Check if empirical standard deviations are coherent with theretical ones
plot(Y_std - Y1_w_std', '.b') % plot empirical differences for Jos function
hold on
plot(Y_std - Y2_w_std', '.r') % plot empirical differences for Gianluca Dorini's function

%%%%%%%%%%%%%%%%%

It looks as though both functions statistically behave very well. They are both excellent functions. Thank you.

Wow, it really is a nice function, very fast and efficient. It is exactly what I was looking for. Thank you.