The Statistics Toolbox implements the gap statistic as a class in the package `clustering.evaluation` since R2013b: http://www.mathworks.com/help/stats/clustering.evaluation.gapevaluationclass.html
load fisheriris;
rng('default'); % For reproducibility
eva = evalclusters(meas,'kmeans','gap','KList',[1:6])
figure;
plot(eva);
You can also use this file exchange: http://www.mathworks.com/matlabcentral/fileexchange/37905-gap-statistics