First of all, the prime geodesic theorem is due to Sarnak, in Margulis' result the geodesics need not be prime. Secondly, for constant curvature, the result goes back to Huber and Selberg, with precise error terms. Margulis extended this to variable negative curvature (with no error term), and that was sharpened (at least theoretically) by Dolgopyat to show that there is an exponentially decaying error term, but the exponent is not explicit at all.
Thirdly, the order of growth result for simple geodesics is due not to Maryam, but to yours truly:
%0 Journal Article
%A Rivin, Igor
%T Simple curves on surfaces
%J Geom. Dedicata
%L MR1866856 (2003c:57018)
Maryam proved an asymptotic result, which (as is not surprising to number theorists) is weaker in practice (my result has explicit bounds, hers does not), but certainly conceptually better (with modern technology - mixing of Teichmuller geodesic flow - one can presumably show the existence of error term a la Dolgopyat, but nothing explicit).
Lastly, and getting back to the actual question, the connection between prime numbers and closed geodesics goes via the Selberg zeta function, and its analogy to the Riemann zeta function. So, I will not speak of prime numbers, but the lattice point counting. It turns out that the counting of all closed geodesics and simple closed geodesics is conceptually very similar, despite completely different result that one gets: counting closed geodesics corresponds to counting points of an orbit of a Fuchsian group (for example, $PSL(2, Z)$) in a disk in the hyperbolic plane (that was what Selberg counted initially). The insight is that simple closed curves lie in a finite set of orbits of the mapping class group, and so counting simple closed geodesic reduces to lattice point counting for the mapping class group on Teichmuller space. I am not sure what the number-theoretic analogue would be (that is, arithmetic geometers certainly study moduli, but I am not sure what easy-to-understand things this counts). Notice that the way Huber/Selberg count lattice points in symmetric spaces is via harmonic analysis, and so what is needed for a parallel theory for Teichmuller space is developing tools for analysis on Teichmuller/moduli spaces. Mirzakhani made the first step in this, by combining Greg McShane's amazing identity (or an extension thereof) with a simple idea of Siegel. This allowed her to integrate the constant function. A baby step, to be sure (an anecdote: I applied for an NSF grant once to develop analysis on Teichmuller space. It was turned down as "too narrow").