function x = mynanvar(x,flag,dim)

xndims = ndims(x);
if dim > xndims
    x(:) = 0;
    return
end
xindex = ~isnan(x);
n = sum(xindex,dim);
x(~xindex) = 0;

% Subtract mean from the data along dim.
xmean = sum(x,dim) ./ n;
x = bsxfun(@minus,x,xmean);

% Compute the sum of squares along dim.
x = sum(x.^2,dim);
if flag == 0
    n = n - 1;
end

% Compute variance whenever number of available data points is
% sufficient.
nindex = n ~= 0;
x(nindex) = x(nindex) ./ n(nindex);
x(~nindex) = 0;

end