You can still use nansum, if you catenate your n-d arrays along the n+1st dimension.
For 2D
% commands de-nested for readability. You can do this with a single line, of course
tmp = cat(3,a,b);
c = nansum(tmp,3);
The general case
function out = nansumByElement(A,B)
%NANSUMBYELEMENT performs an element-wise nansum on the n-D arrays A and B
% A and B have to have the same size
% test input
if nargin < 2 || isempty(A) || isempty(B) || ndims(A)~=ndims(B) || ~all(size(A)==size(B))
error('please pass two non-empty arrays of the same size to nansumByElement')
end
% calculate output
nd = ndims(A); % get number of dimensions
% catenate and sum along n+1st dimension
out = nansum(cat(nd+1,A,B),nd+1);