Bayesian analysis of large geostatistical data faces heavy and expensive matrix calculations. These calculations will be even heavier for multivariate spatial and spatiotemporal data with complex dependency structures. This poses serious problems, such as slow speed and chain convergence, for MCMC sampling algorithms, which are commonly used in Bayesian analysis of spatial models. To avoid such computational problems, an alternative approach is to use low-order models, which improve the convergence rate of MCMC algorithms and the speed of calculations by reducing the parameter space and avoiding heavy matrix calculations. In low-order models, the spatial information of the observed locations is summarized into a set of smaller locations. This smaller set is known as the node set. Determining the node set points and their number such that the estimate of their corresponding spatial dependency structure is a clear and error-free representation of the dependency structure obtained from all the data is a basic and key aspect in building low-order models. Designing the spatial points and number of nodes to implement this dimensionality reduction is the main goal of this paper. To demonstrate the performance of different designs in this class of models, we have analyzed water quality data from a large area of ​​Golestan province during the period 2003 to 2013.