It sounds to me like you are making some confusion between number of variables and order of the integration.
What you are looking for is not a 4th order method, but just a general method to integrate a 4 variables function.
The order of an integrating scheme/method refers to the order of magnitude of the neglected terms (truncation error).
I think matlab does not have any embedded function for integration of functions of more than 3 variables, nonetheless you can write your own integrating algorithm.
Take any numerical analisys book. There you should find integration methods for at least 2 variables. Go up to 4 in a hyper-rectangular domain (that is what matlab handles)is a mere exercise of extrapolation from the latter.
Although traditional methods will work, , they will involve at least 4 for cycles, and this will result in a very slow integration process. What I would suggest you is to check out some Monte Carlo integration schemes; these schemes use a statistical approach to the calculation of the integral and they are intrinsically multidimensional. The only drawback is the low accuracy obtainable in a reasonable time.
