Variable-order fractional advection-diffusion equations have always been employed as a powerful tool in complex anomalous diffusion modeling. The proposed paper is devoted to using the meshless method to solve a general variable-order time fractional advection-diffusion-reaction equation (VO-TF-ADRE) with complex geometries. The proposed method is based on the improved backward substitute method (IBSM) in conjunction with the finite difference technique. For temporal derivative, the finite difference technique and for spatial derivatives, the IBSM are utilized to discretize the equation. The newly developed method is an RBF-based meshless approach, whose solution is constructed by the primary approximation and a series of basis functions. The primary approximation is given to satisfy boundary conditions. Each basis function is the sum of radial basis functions and a specific correcting function. Seven different types of numerical experiments are analyzed to validate the efficiency and wide applicability for multidimensional VO-TF-ADREs.