The main purpose of the natural orbital analysis is to study many particle systems. There, the situation is significantly more complicated because one has several protons and it turns out to be difficult to generalize the decomposition in Eq. 5.15 to the many-proton case because one has no a priori criterion to which proton a particular open path configuration belongs. The problem exists for classical as well as for quantum mechanical protons. However, there is a method that is conceptually correct but not feasible for real applications. For a fixed configuration of protons, one would store the electron two particle density matrix as a matrix where is the number of grid points in a spatial direction. After a sufficiently long MC simulation, one finds a converged results for each matrix element and then can diagonalize the matrix. For low temperature and low density the matrix is block diagonal and the eigenvectors correspond to different localized electronic states, each corresponding to a particular proton. This method shows how one would in principle generate the eigenstates for a many particle system. However, it is not practical because it requires the storage of this enormous matrix and extremely good statistics to fill all relevant matrix elements. The following sections discuss how to approximate this result with computationally feasible methods.