As well known, the positivity is an essential topic when studying gene regulatory networks since the variables involved, e.g., the concentrations of mRNA and proteins, can never be negative. However, the positivity of Caputo fractional order models has been a longstanding problem due to the nonlocality of Caputo fractional derivatives (CFD). In this paper, we present the system comparison method to prove the positivity of Caputo fractional order gene regulatory networks (CFOGRNs) only under positive initial conditions. Moreover, it is found that the positivity results can make it feasible to give proper comparison systems for CFOGRNs, in which the upper and lower estimations can be used to guarantee the stability of the objective CFOGRNs. Thus, the system comparison method for the stability of CFOGRNs is also provided here. Compared to the existing Lyapunov direct method, the proposed system comparison method affords an alternative method for stability analysis and different insights in stability conditions. Finally, these theoretical derivations are illustrated and validated by an example with numerical simulations.