Hi Suzanne,
You're pretty much correct in your first three questions. I'll add a few notes:
1. The values I see in the voxel layer do come from the EBK3D analysis. That the EBK3D didn't show them to me before is because it has made its (and our) life a bit easier (runs quick) by visualising a triangular grid and 'spreading' the range of values of the original data set over it?
Yes. To build a representative histogram of the output surface, the geostatistical layer would have to make many calculations. By borrowing the histogram of the input points, it can process much more quickly. However, it is just a rough visualization of what the surface actually looks like.
2. That means that whatever export format I choose for the statistical layer (for example 'to points') I will obtain the values that I now see in my voxel layer?
Yes, all the geostatistical layer does is call a function that returns a prediction (and possibly standard errors) at a coordinate and then writes the value(s) to an output format. As far as the geostatistical layer is aware, exporting to points, voxels, rasters, etc are all just sets of coordinates used to call this function.
After exporting to whatever format, the renderer of the output is always built on the predicted values stored in that format. So, for example, exporting to a raster will build the render using the entire histogram of the raster. Similarly with a voxel, the renderer is built on the range of the entire 3D volume (with extremes excluded, as Andrew noted).
3. Some of the values I obtain are nonphysical (i.e. the negative values), does this not mean something is fundamentally wrong about the statistical analysis? Especially if the values are really quite different from the original data set (for example if my data ranges from 0.1 to 80 and the voxel layer shows results between -20 and 60)?
Most likely, this behavior indicates that the EBK3D model doesn't fit the data well. Kriging of all types has a tendency to smooth the input data, where the range of the predictions is narrower than the range of the input data. However, it can also make predictions outside of the range of the input points. You're actually seeing both at the same time, where both the min and max of the predictions were lower than the min/max of the original data points; it smoothed off the top but extend the bottom of the distribution. I suspect this is due to the shape of the histogram of the original points, which was heavily right-skewed and had constant values.
You should try to identify where the largest deviations are in the 3D volume. It's possible that you'll find them far off in a corner away from the points, but the still model fits well in the areas where you have data. If that doesn't work out, I don't have very many recommendations aside from what I already suggested: trying different subset sizes, overlap factors, and transformations.
4. Wouldn't adjusting the min/max range in the above mentioned case from -20 and 60 to 0.1 and 80 be cheating? It seems I would then be ignoring what the geostatistical analyst has told me..
I wouldn't call it cheating as long as you disclose what you're doing. It's well understood that kriging doesn't understand what the data means (so it won't respect physical impossibilities), and the predictions become unstable when extrapolating. Simply removing these areas or setting them transparent is very common, but should of course be disclosed. You're not actually changing any data when you do this, just choosing to not display results that you know to be misleading or impossible. Ideally, you'll just find a perfect EBK3D model that doesn't require this, but assuming this doesn't happen, in my opinion, it's better to not display results you know to be wrong.
-Eric
: 
