land = shaperead('landareas','UseGeoCoords',true);
geoshow(land,'FaceColor',[245 237 220]/255)
I want to do this for a 3D globe, but it seems that since 'globe' is not an actual map projection, it doesn't have a frame for which I can set a background color. Thus, when I overlay my data, you can see right through the globe in between the shapes:
h = surf(X*radius,Y*radius,Z*radius,'facecolor',color,'edgecolor','none');
To determine the radius, plot the globe as you have already done, then type xlim or ylim to get the limits of the axes. I don't have Matlab with me right now, but I'm guessing it will be in the range of 0.5 or something. Tinker with the value of the radius until you get a white sphere inside the globe, not too big, and not too small.