The physical interpretation of a variable-order fractional diffusion equation within the framework of the multiple trapping model is presented. This interpretation enables the development of a numerical Monte Carlo algorithm to solve the associated subdiffusion equation. An important feature of the model is variation in energy density of localized states, when the detailed balance condition between localized and mobile particles is satisfied. The variable order anomalous diffusion equations under consideration can be applied to the description of transient subdiffusion in inhomogeneous materials, the order of which depends on the considered spatial and/or time scale. Examples of numerical solutions for different situations are demonstrated. Considering variable-order fractional drift, we calculate and analyze the transient current curves of the time-of-flight method for samples with varying density of localized states.