The study of dual path geometries in the plane goes back to Lie, Tresse, and Cartan. There are some comprehensible modern treatments, but there are also some careless ones, so caution is needed. There is some useful information about dual path geometries near the end of the Ivey-Landsberg book, "Cartan for Beginners".

The basic results are these: If $y'' = \Phi(x,y,y')$ defines a path geometry on a surface $S$ (i.e., a $2$-dimensional family $\Sigma$ of curves such that one curve $\pi\in\Sigma$ passes through a given point $p\in S$ with a given tangent direction), then these curves are the (unparametrized) geodesics of a projective connection on $S$ if and only if $\Phi(x,y,p)$ is at most cubic in $p$, i.e., $\Phi_{pppp} = 0$.

The $3$-dimensional submanifold $Z\subset S\times \Sigma$ consisting of those $(p,\pi)$ such that $p$ lies on $\pi$ is known as the *incidence relation* of the path geometry, and it defines a path geometry on $\Sigma$ by identifying $p\in S$ with the curve consisting of those $\pi \in \Sigma$ that pass through $p$. This is known as the *dual* path geometry. Note that one can identify $Z$ with (an open subset of) the projectivized tangent bundle of either $S$ or $\Sigma$, and that a choice of local coordinates $(x,y)$ on $S$ induces local coordinates $(x,y,p)$ on $Z$ in a natural way.

The condition that the dual path geometry on $\Sigma$ be the (unparametrized) geodesics of some projective connection on $\Sigma$ is that $I(\Phi)=0$, where
$$
I(\Phi) = D^2(\Phi_{pp})-\Phi_p D(\Phi_{pp}) - 4D(\Phi_{yp})
-3\Phi_y\Phi_{pp}+ 4\Phi_p\Phi_{yp} + 6\Phi_{yy}
$$
and where the operator $D$ on any function $G(x,y,p)$ is defined to be
$$
D\bigl(G(x,y,p)\bigr) = G_x(x,y,p) + p\ G_y(x,y,p) + \Phi(x,y,p)\ G_p(x,y,p).
$$

Given local coordinates $(\xi,\eta)$ on $\Sigma$, the dual path geometry is described by some differential equation $\eta'' = \Psi(\xi,\eta,\eta')$, and the condition for it to be the geodesics of a projective connection is that $\Psi(\xi,\eta,\mu)$ satisfy $\Psi_{\mu\mu\mu\mu}=0$ (i.e., $\Psi$ is at most cubic in $\mu$), but $\Psi$ is not directly computable unless you can explicitly integrate the original equation, which is not usually the case. That's why its nice to have the expression $I(\Phi)$, which is explicitly computable directly in terms of $\Phi$.

Interpreting the vanishing of $I(\Phi)$ at a point of $Z$ as 'infinitesimally Desarguesian' may be an attempt to relate this vanishing to something that one can 'see', such as Desargues' configuration. However, I don't think that being 'infinitesimally Desarguesian' at every point implies that the path geometry is actually Desarguesian, so I'm not so sure that this is a good idea, just on terminological grounds. Moreover, I have never seen a precise definition of 'infinitesimally Desarguesian' in the literature, so I don't know where you could find the proof you seek. (However, see below.)

Added Remarks: I thought of a few more things to add that may be of help.

First, given a point $(x_0,y_0,p_0)$ in the domain of $\Phi$,
one can always choose new local coordinates $(x,y)$ in $S$
to bring it to $(x_0,y_0,p_0) = (0,0,0)$ in such a way
that the defining equation becomes $y'' = \bar\Phi(x,y,y')$, where
$$
\bar\Phi(x,y,p) = {\tfrac14}a\ x^2 p^2 + {\tfrac1{24}}b\ p^4 + R_5(x,y,p)
$$
and where the Taylor series of $R_5$ at $(x,y,p)=(0,0,0)$ starts at order $5$. One can't get rid of $a$ and $b$ (if they aren't zero), though one can scale them (they are "relative invariants", so-called). If $ab$ is not zero,
one can't change its sign, which is the same as the sign of
$I(\Phi)\ \Phi_{pppp}$. In fact, the tensor on $Z$ given by
$$
E = I(\Phi)\ \Phi_{pppp}\ (dx\wedge dy\wedge dp)^{\otimes 2}
$$
is a well-defined invariant of the path geometry,
independent of any choice of coordinates,
and its value at the origin in the above normalized coordinate system
is just
$$
E_{(0,0,0)} = ab\ (dx\wedge dy\wedge dp)_{(0,0,0)}^{\otimes 2}.
$$
Thus, there is a distinction between "positively curved"
and "negatively curved" path geometries. I don't know whether one could detect this visually by looking
at some configuration of paths, but it might be interesting to try.

Second, $\Phi_{pppp} \equiv I(\Phi) \equiv 0$ is the necessary and sufficient condition that one be able to choose local coordinates in which the path geometry
is described by $y'' = 0$.

Third, given a point $(p,\pi)\in Z$, you can always choose $p$-centered local coordinates $(x,y)$ on $S$ and $\pi$-centered local coordinates $(\xi,\eta)$ on $\Sigma$
so that the incidence relation of the path geometry defining $Z\subset S\times \Sigma$ takes the form
$$
y = \eta + \xi\ x + {\tfrac1{48}}b\ \xi^4x^2 + {\tfrac1{48}}a\ \xi^2x^4
+ R_7(x,y,\xi,\eta)
$$
where the Taylor series of $R_7$ at $(x,y,\xi,\eta)=(0,0,0,0)$ has no terms of order less than $7$. This normal form of the relation highlights the symmetry between the dual path geometries. (The appearance of $a$ and $b$ in this formula and in the formula above is not a coincidence, nor is the fact that the $a$ and $b$ appear in reverse order to that in the formula above.)

In fact, this formula makes me wonder about Arnold's exercise. I have heard someone who should know (Simon Gindikin, in fact!) say that 'infinitesimally Desarguesian'
means that "the paths agree with straight lines up to third order". The above formula indicates that this happens when $b=0$, i.e., when $\Phi_{pppp}$ vanishes, not when $a=0$,
which would be when $I(\Phi)$ vanishes. I'm afraid that this mystery will only be cleared up when we get a precise definition of 'infinitesimally Desarguesian'. (See below.)

Fourth (and finally), though it may be immodest of me to point this out, there is a treatment of dual path geometries and their corresponding ODEs in my paper with Phillip Griffiths and Lucas Hsu, "Toward a geometry of differential equations", in Geometry, topology, and physics for Raoul Bott, Conference Proceedings and Lecture Notes in Geometry and Topology Volume VI, International Press, Cambridge, MA, 1995, pp. 1–76. (MR 97b:58005) We discuss this in a couple of places, first in Section 1 and then in the Appendix, where we give an exposition of the equivalence problem applied to dual path geometries.

More Added Remarks:

I have thought more about what *infinitesimally Desarguesian* might mean, and I have a proposal that might work:

Say that a path geometry on a surface is *infinitesimally Desarguesian* if, for any path $\pi$ of the family and any point $p$ lying on $\pi$, there exist $p$-centered coordinates $(x,y)$ on the surface such that $\pi$ is defined in these coordinates by $y=0$ and such that the (local) defining differential equation $y'' = \Phi(x,y,y')$ of the path geometry satisfies the condition that $G(x,0,0)\equiv0$ for all of the $G(x,y,p)$ on the following list:
$$
\Phi, \Phi_y, \Phi_p, \Phi_{yy}, \Phi_{yp}, \Phi_{pp}
$$
(i.e., the Taylor series of $\Phi$ with respect to $y$ and $p$ starts at order $3$).
Note that, by the formula for $I(\Phi)$, this hypothesis is sufficient to conclude that the function $G = I(\Phi)$ also satisfies $G(x,0,0)\equiv 0$. In particular, since the vanishing of $I(\Phi)$ is independent of coordinates, it follows that, if the path geometry is 'infinitesimally Desarguesian' at *each* pair $(p,\pi)$ then $I(\Phi)$ vanishes identically. Consequently, the dual path geometry is projective, just as was claimed in the sources cited by Jan.

Note that it's *not* enough to assume that the path geometry is 'infinitesimally Desarguesian' at one pair $(p,\pi)$. One needs to assume it everywhere. Also, note that the condition is *not* that all these partial differential expressions should vanish in *any* local coordinate system on the surface, which is impossible, but that, for each pair $(p,\pi)$, there exists *some* coordinate system in which these expressions all vanish.

Now, the converse, which would be that, if $I(\Phi)\equiv0$, then for every incident pair $(p,\pi)$ there exists a $p$-centered coordinate system $(x,y)$ in which $\pi$ is defined by $y=0$ and for which $G(x,0,0)\equiv0$ for all the $G(x,y,p)$ on the list
$$
\Phi, \Phi_y, \Phi_p, \Phi_{yy}, \Phi_{yp}, \Phi_{pp}\ ,
$$
seems to be a bit messy to verify, but does appear to be true. I just don't see why anyone would want to make this definition, when the 'real' definition should be that $I(\Phi)$ vanishes.

Finally, it's not as though the paths really closely approximate straight lines in such coordinates. One *doesn't* get that the paths look like straight lines to 3rd order in such coordinates. (However, the dual path geometry will have this property.) I think that, if $\Phi_{pppp}(0,0,0)$ isn't zero, it's not possible to make the paths near the corresponding $(p,\pi)$ all look like straight lines up to 3rd order near that point. What you *can* say is that, in these coordinates, one has the incidence relation in the form
$$
y = \eta + \xi x + P(\xi,\eta)\ x^2 + O(x^3)
$$
where $P$ vanishes to order at least $3$ in $(\xi,\eta)$.

Thus, the terminology and motivation for this (conjectural) definition of *infinitesimally Desarguesian* is something of a puzzle, to me. I guess the proposers wanted something to state that would capture $I(\Phi)\equiv 0$ without having to actually define $I(\Phi)$.