Systematic sampling is one of the simplest, easiest and the most common used sampling methods. However, when the population size N is not a multiple of the sample size n, the systematic sampling cannot be performed. Not only is it difficult to determine the sampling interval k, but the sample mean will be a biased estimator of the population mean. To solve this problem, this paper proposes an improved method for the systematic sampling: the remainder systematic Markov chain design. The first- and second-order inclusion probabilities are derived, yielding the Horvitz-Thompson estimator and its variance. The simulation results demonstrate the effectiveness of the proposed method for different super-populations.