It's definitely possible, though it doesn't have a clean explicit formula. The projection of a sphere is an ellipse, and given the projection matrix you should be able to find an explicit formula for the ellipse, something of the form ax^2+bxy+cy^2+dx+ey+f=0 (the canonical formula for a conic section); from this you can find a rational parametrization of the ellipse (of the form x=f(t), y=g(t) with f and g rational functions). The problem then comes down to minimizing (f(t)-x_0)^2+(g(t)-y_0)^2 as a function of t, and it turns out that rootfinding methods work pretty well for this - there are exact solutions but they require solving a degree-4 equation, which can involve hairy complex arithmetic and tends to be less numerically stable (and thus less accurate) than the approximate methods anyway.