You need to compute the total population within all possible distances, plot them together, and interpolate within that plot.
I will illustrate using population for the District of Columbia (USA) represented on a 10m grid. It begins by computing the population density (the units do not matter; I used people per acre) and converting that to a grid. For the purposes of the density computation, obtain accurate areas (perhaps by projecting the original data using an equal-area projection). For the conversion to grid format, however, use a projection suitable for accurate calculation of distances.
More intense colors denote higher population densities.
Next, compute the Euclidean distance grid to your points. Here is a map of the population density grid hillshaded according to the distance grid to a few hundred randomly located points. (You cannot see the points at this resolution, but their locations are evident as the cusps of the "valleys" in the hillshade.)
Population densities graduate from blue through green to yellow.
Third, discretize the distance to create zones to summarize. This procedure is done by dividing the distance by the zone width and truncating (using the
floor functions). A zone width less than the cellsize is unnecessary, and likely a much greater width will yield acceptable accuracy and improved calculation times, so start with a large width--maybe a few kilometers at a country scale--and narrow that if the results are too imprecise. For this example I used a width of 10 meters (the cellsize). (No illustration is needed because the zone grid looks almost exactly like the distance grid.)
Using these distance zones, compute a zonal sum (or mean) of the population density grid. The output will be a table containing at least two columns: zone identifier and total population density within that zone.
The total population density is directly proportional to total population. Using a spreadsheet or statistical software, then, compute the cumulative sum of the zonal sum (ordered by distance) and rescale that to the range 0 - 100% (just divide the cumulative sums by the total sum). The plot of distance against the proportion of population provides the desired information.
In this table,
Value is the zone identifier: 0 for 0-10 meters, 1 for 10-20 meters, and so on. This was multiplied by 10 and then 5 was added to produce the calculated
Distance column representing a typical distance within each zone.
For example, the tenth percentile of distance is found by reading the horizontal axis to 10% and looking up to the graph, which is at a height between 125 (9.45%, on line 14) and 135 meters (11.02%, on line 15). For greater precision, linearly interpolate between these values: the target of 10% is 10 - 9.45 = 0.55% greater than 9.45% and 11.02 - 10 = 1.02% less than the two surrounding values, which are 11.02 - 9.45 = 1.57% apart. The interpolated distance therefore must be both 0.55 / 1.57 times 10 meters greater than 125 and 1.02 / 1.57 times 10 meters less than 135; both calculations give 128 meters. (Any values after the decimal place would be meaningless, given that the distance of each zone already represents a 10 meter range.)
Similarly, we can read off the other percentiles:
Read the rows like this: "10% of all people in this dataset are located within 128 meters of a point; 25% are located within 215 meters of a point;..." and so on.
Comments and advice
For working at a country wide scale, these calculations can be tiled into smaller grids. Carry out the work for each tile (using distances to all points, whether or not they are located within the tile) up to the point where the zonal summary table is obtained. Afterwards, combine all the zonal summary tables by adding their values across rows: for each distance zone, there will be one sum for each tile and adding these gives the total population in that zone for all tiles. Generate the percentiles from this combined table.
One advantage of the tiled approach is that a different projection can be used for each tile, allowing more accurate calculation of distances. To do this, create the tiles using a fixed projection for the country, then reproject each tile (and all the points) to its local coordinate system for doing the distance calculations and zonal summaries. You may safely throw away these reprojected grids after they are used.