Rectification of the ancient geographic coordinates in Ptolemy ’ s Geographike Hyphegesis

Abstract. A multitude of the ancient places given by Ptolemy in his Geography (~150 AD) are so far unknown. One of the main problems of their identification are the errors of the given ancient coordinates. The different kinds of errors are illustrated by examples. A new geodetic-statistical analysis method is described, by which groups of places with homogeneous systematic errors and places with gross errors can be determined. Based on a transformation function describing the systematic errors, presumable modern coordinates of unknown places can be computed. That, in conjunction with further information, can make possible their identification. A test of the analysis method is carried out on a complex simulated example and shows its practicability. The analysis method has been applied within an interdisciplinary research project on Ptolemy's Geography. Further developments of the method are imaginable to make it accessible for related data diagnostics.


Introduction
The development of information technology has offered new possibilities for the history of geography and cartography, which have found their way into scientific investigations in the last decades (e.g.Beineke, 2001;Niederöst, 2005;Kleineberg et al., 2010).In this context, the fields of interest are especially the methods of historical measurement and cartography including their accuracy and defects as well as the interpretation of historical information specifying the position of unknown historical places.
In the present work an analysis method is described, which has been developed and applied to investigate the geographic coordinates in Claudius Ptolemaios' (Ptolemy) Geographike Hyphegesis (Geography, ∼150 AD), one of the most important ancient scientific works.In addition to a theoretical part on the cartographic methods of mapping the world, the Geography contains a catalogue of locations (Books II-VII) with geographic coordinates of over 6300 ancient places (towns, settlements, rivers, mountains, border points, and others) covering the Ecumene, the then inhabited world known to the Greeks and Romans.A multitude of them have not been identified to date.Key to the modern counterparts are the an-cient coordinates; however, due to their errors, their direct usage is not possible.Further, the question of the identification is hindered by the fact that the original Geography does not exist any more; it has been handed down by Greek medieval manuscripts, which show differences.The manuscripts are presumably based on two recensions, Ω and Ξ, which were revised and published by Stückelberger and Graßhoff (2006).
So far there is a lack of an appropriate computational method, by which the errors (distortions) of the ancient positions can be analyzed and rectified.Rare efforts in this direction are: For regions of Books II and III Cuntz (1923) compares the distances between the Ptolemaic places, derived from their coordinates, with distances of itineraries.Strang (1997) explains the distortions of Hibernia (II.2) and Albion (II.3) by scaling errors and partly by rotations.Tsorlini and Livieratos (2007) and Tsorlini (2009) work on Hispania (III.4-6),Southeast Europe , and Asia Minor (V.2) by computer based visualization methods for error diagnostics and use for a rectification a polynomial transformation function.
In recent years an interdisciplinary project group at the TU Berlin has worked on the Ptolemaic data with the aim to identify unknown Ptolemaic places on the basis of their ancient Published by Copernicus Publications.
C. Marx: Rectification of the ancient geographic coordinates in Ptolemy's Geographike Hyphegesis coordinates.Due to the errors and inaccuracy of the coordinates, further information is consulted for the identification, that is topographic and historical information, archaeological sites, other ancient texts, and toponymy.The findings of the project are not the last word on the subject; in fact, specialists such as historians and archaeologists are asked to confirm or disprove them.
Concerning the mathematical questions of the project, the problems and aims are to determine systematic distortions, gross errors (scribal errors of the manuscripts), and the local accuracy of the ancient data as well as incorrect modern identifications and modern coordinates for unidentified ancient places.For this purpose a geodetic-statistical analysis method was developed and implemented in software.The present work gives a detailed representation of the actual analysis method, which is described only briefly in Kleineberg et al. (2010Kleineberg et al. ( , 2012) so far.Further, the errors of the Ptolemaic coordinates are illustrated and a simulated example is given, by which the functioning of the analysis method is demonstrated.

The errors of the Ptolemaic coordinates
For the places in the catalogue of locations Ptolemy gives the geographic longitude Λ and latitude Φ.The difference between the Ptolemaic and the modern geographic coordinate system is the zero meridian.Ptolemy locates the zero meridian at the "Blest Islands" (IV.6.34).The modern counterpart of that isles is disputed; two possibilities come into consideration, the Canaries (e.g.Tsorlini and Livieratos, 2007;Stückelberger and Mittenhuber, 2009, p. 240) as well as the Cape Verde Islands (e.g.Reichert, 2003;Rawlins, 2008, p. 573).Regardless of which zero meridian is taken into account, there remain large deviations between the converted ancient and the modern longitudes.Also, ancient and modern latitudes often differ considerably.Regarding the identification of places, it is of interest whether the coordinate differences contain systematic components.If so, they can be described by a mathematical function, and the positions of unidentified places can be rectified.
To illustrate the occurring deviations between Ptolemaic and modern coordinates, exemplarily 46 places in Illyricum (II.16; situated between the Adriatic and the Dinarides) are considered, which are certainly identified (see Kleineberg et al., 2012, pp. 189).The following computations are based on the Ω-recension.The differences between modern and Ptolemaic longitudes amount to 23 • to 29 • .Thus, for this example the Cape Verde Islands are assumed to be the location of the zero meridian, and the ancient longitudes were converted by Λ = Λ + λ SA , where λ SA = −25 • 10 is the modern longitude of Santa Antão, the most western island of the Cape Verde Islands.The differences ∆λ = λ − Λ and ∆φ = φ − Φ with respect to the modern longitude λ and latitude φ were computed.The difference vectors (∆λ ∆φ) are shown in Fig. 1.They have varying directions and often magnitudes of 2 • to 3 • in longitude and 1

2
• to 1 • in latitude.However, there are also vectors which are similar to each other.

That indicates systematic errors (distortions).
A better fit between ancient and modern coordinates can be obtained by a transformation of coordinates, which absorbs the distortions of the ancient coordinates.In general, the transformation function is Λ = f (u,λ,φ), Φ = g (u,λ,φ), where u is the vector of the distortion parameters u j .Using the given places of Illyricum, a computation of u was performed by an adjustment (see Sect. 3) exemplarily for a scaled and shifted transformation where m λ and m φ are scales and Λ 0 and Φ 0 shifts.The parameters of the inverse transformation with respect to Eq. ( 1) are obtained by the inversions Then, the approximately rectified ancient coordinates are The remaining errors of this transformation are the residuals ∆λ = λ − λ and ∆φ = φ − φ.The adjustment of u was carried out using all 46 places.The residual vectors (∆λ ∆φ) are plotted in Fig. 2. They have magnitudes comparable to the vectors in Fig. 1, because the determined transformation function is not well fitting.In Fig. 2 similar residual vectors in local groups reveal shifts of groups of places against each other.The cause of shifts can be reference places, which the localization (by measurement or calculation) of neighbouring places was based on.(Concerning Ptolemy's usage of itineraries, also Cuntz (1923, pp. 112) assumes reference places in Gallia, Italia, and Pannonia Inferior.)When the relative positions of the reference places were inaccurate, the errors were transferred to the associated places.
In addition to the shifts, scalings can be assumed, which appear more wide-ranging.Reasons for scalings can be Ptolemy's overestimation of the longitudinal dimension of the Ecumene (180 • to the most eastern place in China instead of about 130 • ) as well as differences between ancient measurement units (different types of stades), which were unintentionally not considered.(Due to the interactions of different influences, the scalings of the longitudes and latitudes are possibly not entirely identical.) The systematic scalings and shifts can be described by Eq. (1).Groups of places with homogenous scales and shifts (transformation units) and their distortion parameters can be determined by means of geodetic-statistical methods (see below).Subsequently, the coordinates of so far unidentified ancient places can be rectified as an important step towards their identification.
Residual vectors in Fig. 2, which are similar to each other, have an apparently random component in addition to the common systematic component.Data used by Ptolemy, that is measurement data, data derived from measurements, and information from travel reports (itineraries), necessarily had random errors, which propagated to the final coordinates.Also the combination of data by Ptolemy led to random adulterations.For example, Ptolemy reduced given distances along routes, which were not straight but curvilinear, by a unitary factor (Berggren and Jones, 2000, p. 17), whereby the real distances were enlarged or shortened.• .A study on the frequencies of the fractions showed that the coordinates are given with different resolutions ≥ 5 and that there are regional differences in this respect (Marx, 2011).Due to the rounding of coordinate values, random adulterations were induced.
For Illyricum Fig. 5 shows the determined transformation units as well as the residual vectors after a transformation with individual shifts for transformation units.Apart from a few exceptions, the residuals are considerably smaller than in Fig. 2. Their magnitudes are reconcilable with the assumed local accuracy of the Ptolemaic coordinates.Table 1 gives the average of the estimated standard deviations s Λ and s Φ of the ancient longitudes and latitudes (scale-corrected).Their order of magnitude is about 10 km to 20 km (in Il11 s Λ is underestimated due to a low number of places; in the case of a possible underestimation, an average a priori standard deviations are given in Kleineberg et al., 2010Kleineberg et al., , 2012)).
Partly the transformation units in Fig. 5 overlap somewhat.That also occurs in other regions of the Geography.Possibly, in these cases the positions of neighbouring places were determined independently.Coordinate values which cannot be explained by systematic scalings and shifts are regarded as grossly erroneous.Reasons for such errors are the following.The copying of the manuscripts caused scribal errors and intensional changes of coordinates are supposable, which may have led to corrections or significant deteriorations.Further, possibly the data underlying Ptolemy's work were erroneous, Ptolemy made mistakes, or he only gave crude positions in case his information was incomplete.
Based on assumptions about the ancient coordinate accuracy, gross errors can be detected by statistical hypothesis tests.The problem is aggravated by the fact that, dependent on the region, more or less places have no certain modern identification but more than one suggestion for the identification.Regarding the distortion model in Eq. ( 1), a wrong identification farther away from the correct place can be considered as a gross 2-dimensional model error.In Fig. 5 there are five places with residual vectors differing from their neighbouring vectors, which show gross coordinate errors.Two of them are wrongly located islands, two errors can be explained by another ancient coordinate variant and one error can be explained by a scribal error (Kleineberg et al., 2012, p. 191).
Due to the interaction of the various errors and the varying accuracy of the coordinate values, it is possible that the relative position of two places is contradictory (e.g.too large a distance or western instead eastern).Nonetheless, their coordinates can be acceptable in case there are other places with the same systematic errors.An example of contradictory ancient positions are the strong distortions in Gallia recognizable by the considerable differences of neighbouring residual vectors (cf.Fig. 3).

The method of rectification
If the systematic errors of the Ptolemaic locations with existing identifications are determined, it will be possible to compute presumable modern coordinates of unidentified Ptolemaic locations.In doing so, also uncertain identifications can be supported or rejected.Further, after a deduction of the systematic errors, the remaining errors reveal the local coordinate accuracy.
In the last decade several efforts have been made to investigate the accuracy of historical maps or to rectify them (e.g.Beineke, 2001;Niederöst, 2005).Common methods applied are the adjustment of parameters of a transformation between coordinates of the historical map and modern coordinates as well as, based on this, georeferencing, a geometric rectification as it is used in remote sensing for the rectification of satellite images (e.g.Albertz, 1999, pp. 93).These methods are not suitable for a rectification of the Ptolemaic coordinates since the problems differ significantly, as it is pointed out in the following.
In the case of georeferencing of historical maps, the digitized map is rectified on the basis of (unique) control points, whose coordinates were determined in the map and the modern reference system beforehand.The parameters of a (continuous) transformation function are determined by an adjustment (one-step procedure); then the rectifying transformation is applied.Generally, the chosen transformation function (mainly similar, affin, polynomial, or rubber sheeting) does not specify the (unknown) causations of the complex distortions and is applied to the whole map.The adjustment is not based on a stochastic model (accuracy of the determined coordinates).Gross erroneous areas of the map are only identified visually.If the transformation function is not realistic, areas with gross and systematic errors (distortions) cannot be identified reliably.The map is not segmented in such areas.The common measures of the accuracy derived from the remaining errors after the adjustment depend on the chosen transformation function and are affected by the existing but disregarded systematic and gross distortions.Then these measures do not represent the accuracy of the sources underlying the map (e.g.historical measurements).As an extension of the typical approach, for historical maps Beineke (2007) suggests a calculative consideration of the remaining errors by a multiquadratic interpolation by Hardy (1972).
In the case of the data of the Geography, the coordinates of single places are given, whose modern coordinates are known (certainly identified), unknown (not identified), or uncertain and partly not unique (> 1 identifications).Places with certain or uncertain identifications can be used as control points for an adjustment.In order to detect wrong identifications and gross errors, the transformation function has to specify the real distortions and a stochastic model has to be used, describing the assumed accuracy of the coordinates.Since the ancient positions have local shifts in groups, the transformation function cannot be continuous (regarding more than one shifted group).Places with homogenous systematic errors and gross errors have to be identified.Finally the accuracy of the ancient coordinates is of interest.It can be assumed that the ancient places have remaining errors of a few arc minutes up to ≈ 1 2 • .The result of the rectification is satisfying if the remaining errors distribute randomly by the majority and have magnitudes explicable by the local coordinate accuracy.An interpolation of the remaining errors would be inappropriate (in terms of a rectifying transformation), because random errors must not be transferred unreservedly to the neighbouring places.Due to the determination of new identifications, the analysis is a repeated process.
In order to detect gross model and data errors for given data, resistent or robust estimations are usual in data analysis (e.g.L 1 -norm, Huber estimator).However, these estimators are not reliable in the presence of a large portion of gross errors.Then a promising approach is a combinatorial analysis of the observations.Here, such a method is used, which Neitzel (2004) suggests for the deformation analysis of geodetic networks (method of the maximum subsample or MSS-method).Depending on the given data, this method can be impractically due to its high computing time.Therefore, the method is extended here by further analysis steps.Moreover, it is adapted to the peculiarities of the present problem.
In general terms, the analysis method developed is a multistage combinatorial search for sets of Ptolemaic places with homogenous distortion.Its procedure is: 1. Initial solution: analysis of the precision of the coordinate values, determining approximative values for the scales, generating initial subsets of places with similar distortions 2. Modified MSS-method: searching for consistent transformation units in the initial subsets 3. Forward-strategy: searching for best possible mergings of unassigned places with transformation units 4. Verification of the scales: testing the suppositional scales introduced in step 1 for validity by an adjustment of scales 5. Merging of transformation units: testing neighbouring transformation units for whether it is possible to merge them 6.Postprocessing: visually checking and manually improving of the results The generation of initial subsets of places in step 1 is necessary because of the computationally intensive modified MSS-method in step 2 and benefits a reliable partition of the places into groups of homogenous distortion.
In steps 1 and 3 more than one ancient coordinate variant per place can be processed.By default the main recensions Ω and Ξ published by Stückelberger and Graßhoff (2006) as well as the editions of the Geography by Nobbe (1843Nobbe ( -1845) ) and Müller (1883Müller ( /1901) ) are used.
Since the random errors of the ancient coordinates may be large, the adjustment of a scaled and shifted transformation, Eq. ( 4), may result in implausible parameter values.Also disadvantageous geometric configurations may cause unrealistic results.Thus, restrictions in form of constant, hypothetical values for m λ and m φ are introduced from step 2. This has also the advantage that the computing time is reduced massively in step 2 (see Sect. 3.3).A test of the hypothetical values is performed in step 4. If the test results in a rejection of a hypothetical value, the analysis has to be repeated with a modified value.
The analysis method is applied to separated investigation areas for which uniform scales are assumed.The areas correspond to the chapters of the Geography or are summarizations of chapters in case correlated distortions are to be expected due to a close position (e.g.Hispania in II.4-6).
www.hist-geo-space-sci.net/3/99/2012/ Hist.Geo Space Sci., 3, 99-112, 2012 If there remain places inconsistent in the determined transformation units, reasons and actions can be: (1) The place is consistent in a neighbouring investigation area; it will be tested there.(2) The coordinate accuracy is lower than presumed; if justified, the a priori accuracy will be decreased.
(3) The identification is wrong; a transformation, Eq. (3), will be performed.(4) A coordinate is grossly erroneous; a scribal error will be quested.If input parameters are modified or identifications are changed or join, the analysis procedure has to be repeated (not necessarily all steps).
In the following the adjustment model and statistical tests applied are described as well as the steps of the analysis procedure.

Adjustment model and statistical tests of consistency
The distortion parameters u j are determined by least squares adjustment (see, e.g.Niemeier, 2002;Jäger et al., 2005), whereby normal distribution is assumed for the random errors.The assumptions about the distortions in Sect. 2 lead to the following model of an adjustment of observation equations.The Ptolemaic coordinates (Λ, Φ) are erroneous, observed quantities, the modern coordinates (λ, φ) are errorfree constant quantities.The functional model (observation equations) of p locations is where i = 1... p and v Λi , v Φi are the corrections (here also referred to as residuals in contrast to textbooks, where often residual = −v).The four unknown parameters m λ , m φ , Λ 0 , Φ 0 contain local and global effects.Especially, m λ contains Ptolemy's underestimation of the circumference of the earth.Λ 0 contains the difference between the Ptolemaic and the modern zero meridian.(Since the Ptolemaic zero point is not known for sure, the difference is not taken into account beforehand; however, that is irrelevant for the results.)m λ and m φ are set constant in steps 2, 3, 5, and 6 of the analysis method.
The stochastic part of the model describes the accuracy and the correlation of the observations by means of the covariance matrix where the diagonal elements σ 2 i of C are the variances, the off-diagonal elements σ i j are the covariances, σ 2 0 is the a priori variance of unit weight, and P is the weight matrix.There is no information about correlations between the coordinates of places as well as between Λ and Φ of a place so that no correlations are applied (C ΛΦ and P are diagonal).
Since the Λ i and Φ i have no common parameters and correlations in the adjustment model, the system in Eq. (4) consists of two independent systems for the Λ i and Φ i , respectively.Assumptions about the accuracies of Λ i and Φ i can be derived from the following considerations.The accuracy is related to the resolution a of the coordinate values.Measurings of the shadow length by a reproduction of a gnomon resulted in standard deviations of the derived latitudes of 2.4 ( Lelgemann et al., 2005).Based on this, a somewhat larger value, such as σ = 4 , is a reasonable value for the standard deviation of latitudes with the most accurate resolution of a = 5 .This is also assumed for the longitudes here.From coordinates rounded to full degree a maximal random error of about 2 3 • can be expected.Using the 3σ rule for gross errors here, according to which the maximal accepted correction v is 3σ, a maximal error of • has a = 5 or a = 10 ); however, at least the portions of different resolutions can be estimated according to Marx (2011).The analysis method is started with small σ i ; after a run only individual σ i are increased if necessary and justified by σ(a).
The solution for the unknown parameters and other quantities is derived by the known formulas, see, e.g.Niemeier (2002, p. 117).
As usual, it is tested whether the adjustment model (functional and stochastic) is valid for a group of observations by the model/global test (e.g.Niemeier, 2002, p. 147).Furthermore, single observations Λ i and Φ i as well as the position vectors (Λ i Φ i ) are tested for gross model or data errors by means of single tests.If a local group of places satisfies the tests, it is considered as homogenously distorted (scaled and shifted).
Test statistic of the model test is (e.g.Niemeier, 2002, p. 150) where r is the redundance, S is the weighted sum of squared residuals, s 2 0 is the a posteriori variance of unit weight, and v is the vector of the corrections v i .T G will be χ 2 r -distributed, if the adjustment model is correct (one-sided test).Test statistic of the one-dimensional single test with uncorrelated observations is the standardized residual (e.g.Jäger et al., 2005, p. 193) where p i is the weight, r i is the partial redundance of the observation, and s vi is the standard deviation of the correction v i from the adjustment.w i will be N(0,1)-distributed, if there is no model error with regard to the i-th observation (two-sided test).From a general test statistic for a group of observations (Jäger et al., 2005, p. 189) a test statistic for a position (Λ i Φ i ) can be derived.Due to the independence of the two adjustment systems for the Λ i and Φ i , the result is T P will be F 2,∞ -distributed (and 2T P χ 2 2 -distributed), if there is no model error with regard to (Λ i Φ i ) .
By his B-method Baarda (1968) adapts the significance level α G of the model test to the significance level α S of the single test to make both tests comparably sensitive.To this end, the powers of both tests are set to β G = β S and the noncentrality parameters to λ G = λ S .Since the model test can become too sensitive, a modification by Hahn et al. (1989Hahn et al. ( , 1991) ) is used here, where it is set λ G = kλ S with 1 < k < 2.
The size of the error of an observation can be estimated by (e.g.Jäger et al., 2005, p. 194).

Coordinate accuracy
In Marx (2011) a method is given to estimate the proportions of different coordinate resolutions for a considered group of places.The method is applied here; however, in an investigation area with a lower numbers of places the results give only an approximate insight into the occurring precisions.According to the results, initially a unitary standard deviation is chosen by the model in Table 2.

Scales
Places of the investigation area with certain or uncertain identification(s) are used as control points for a first adjustment of average parameters of distortion model Eq. ( 4).Since the shifts within the investigation area are not considered, the resulting scales will be adulterated.However, the influences of the different shifts will diminish each other so that on average the entire influence is only low.If there are reasonable ancient coordinate variants and/or identification variants, the consideration of all variants can lead to a better estimation result in comparison to a possibly choice of wrong variants.In doing so, for each place all combinations of ancient and modern variants are generated.Each combination presents an observation equation and its influence is weighted down so that the influence of all combinations of a place equals the influence of a place without variants.
A further way to gain more reliable results is the use of the initial subsets of places (see Sect. 3.2.3),whereby the existing shifts are considered approximately.Then, in the adjustment model there are the joint scale parameters as well as N groups of places with their own shift parameters.The observations equations are where k = 1...N, i = 1...n k , and n k is the number of places in the k-th group.

Initial subsets of places
Based on the results of the adjustment of average parameters of distortion a residual plot is generated, showing v Λi and v Φi .Although the adjusted scales are only approximate values, groups with homogenous shifts are indicated by similar vectors (v Λi v Φi ) .Consequently, the places are grouped into initial subsets of homogenous shifts.In the case of ancient coordinate variants and/or identification variants one variant is chosen (Ω by default).

Modified MSS-method (step 2)
The aim of the MSS-method by Neitzel (2004, p. 109) is to find the maximal subset of the given data whose least squares adjustment has an agreeable result.The search is performed by going through all reasonable combinations of the observations.Combinations whose adjustment result is not agreeable are discarded.The method is adapted to the present problem as follows.
As criterion for the agreeability (consistency) the model test Eq.( 6) and the single tests Eqs. ( 7), ( 8) are used here.With the increasing portion of inconsistent data the MSSmethod becomes more computationally intensive and finally impracticable.Thus, the method is applied only in the initial subsets of places so that the data volume is relatively small.Furthermore, the number of combinations can be reduced massively here by a limitation to reasonable ones (see below).In the case of more than one maximal consistent subset the weighted sum of squared residuals S is used here as a criterion.
To make the method applicable in areas with a large portion of uncertain identifications (often > 1 alternatives), it was modified such that more than one identification per place can be processed.For it, the different identifications of an ancient place are considered in the combinatorial run.That can increase the computational time enormously so that ancient www.hist-geo-space-sci.net/3/99/2012/ Hist.Geo Space Sci., 3, 99-112, 2012 coordinate variants are not considered in this step.From experience, the Ω-recension often provides a proper variant so that this variant is used by default.Since the adjustment result could degenerate with variable scales m λ and m φ , both parameters are treated as constants and get approximative values (see Sect. 3.2).Minimum constraint of model Eq. ( 4) is (since (4) consists of two independent systems without correlated observations, cf.Sect.3.1).Due to the constant scales, the observation equations ( 4) can be rearranged to By means of the substitutions l Λi = Λ i − m λ λ i and l Φi = Φ i − m φ φ i the observation equations of the mean are obtained, whereat l Λi and l Φi are the observations of the mean parameters Λ 0 and Φ 0 .The minimum constraint of model Eq. ( 12) is also (11) so that Eq. ( 12) is equivalent to Eq. ( 4).Applying the MSS-method to determine the mean, the number of combinations of the observations can be reduced by sorting the observations by size; only successionally sorted observations compose a reasonable subset (Neitzel, 2004, p. 110).At first, the search for consistent subsets is performed for longitude and latitude separately.Tests applied are the model test Eq.( 6) and the single test Eq.( 7).From the consistent subsets found for longitude and latitude the best maximal consistent subset of locations is derived.For that the combinations of the subsets in longitude and latitude are gone through and for a combination the 2-dimensional single test Eq.( 8) is applied to the involved locations.The best maximal subset of locations found forms a transformation unit.Among the remaining locations the search is continued.

Forward-strategy (step 3)
As a result of the application of the MSS-method in the initial subsets of places there are N transformation units U k (k = 1...N) and in general further places without association to a unit.Possibly, among them there are places with correct coordinates and identifications.
In deformation analysis of geodetic networks the typical problem is to decide, based on a group of stable points, whether further points are deformed or not.Then a forward strategy is common (e.g.Niemeier, 2002, p. 380), where the questionable points are involved sequentially in the congruence hypothesis (statistical test).A similar approach is used here.
Causes of a remaining unassigned place can be: (1) consistency in a transformation unit of another initial subset of places, (2) grossly erroneous ancient coordinate variant(s) used in the MSS-method, (3) wrong identification(s), (4) accuracy in stochastic model too high.Causes (1) and ( 2) are tested, that is, whether an unassigned place is consistent in a transformation unit with one of its identifications and ancient coordinate variants.
For a given identification the possible transformation units are limited to the neighbouring ones (modern system).Two tests are applied: 1.The convex hull (generated by Jarvis march, e.g.Preparata and Shamos, 1985, p. 110) of the considered transformation unit is enlarged by a buffer zone with given width b so that a polygon is generated which represents the unit and its neighbourhood.By means of a point-in-polygon test (carried out by the Jordan curve theorem, e.g.Bill, 1996, p. 28) is tested, whether the identification lies within the polygon.
2. It is tested, whether the spherical distance between the center of the transformation unit and the identification is lower than a given maximal value d.
For each unassigned place the given identifications are gone through and the neighbouring transformation units are determined.For each pair of identification and transformation unit the ancient coordinate variants of the place are gone through.The considered unit is enlarged by the place and adjusted, then the locations are tested by the single tests Eqs. ( 7) and ( 8).Among all possible mergings of places and transformation units the merging with the lowest 2-dimensional test statistic Eq. ( 8) is selected.The process is repeated with the remaining unassigned places.
Finally the model test Eq.( 6) is applied on all changed transformation units.

Verification of the scales (step 4)
In the adjustment model of the MSS-method (step 2) and the forward-strategy (step 3) the scale parameters are set to hypothetical values and are constants.On the basis of the formed transformation units U k (k = 1...N), the scales can be computed much more reliably than by a small number of places.By a t-test is tested, whether the computed scales equal the hypothetical values.
The adjustment model corresponds to Eq. ( 10), in which there are the joint scale parameters as well as N groups of places U k with their own shift parameters.Test statistic of the t-test (e.g.Niemeier, 2002, p. 66, 356) for one scale parameter is where m (0) is the hypothetical value, m is the adjusted value, and s m is its standard deviation.T m will be t r -distributed, if the expectation value E(m) of m is the hypothetical value m (0) (two-sided test, r is the redundance).
If there is a significant difference between an adjusted and a hypothetical scale, the hypothetical scale has to be modified and the analysis method has to be repeated.

Merging of transformation units (step 5)
By the generation of initial subsets of places with homogeneous distortions (step 1) groups of places are possibly separated which have no significantly different distortions.It is necessary to merge such groups.To this end, a search for maximal (in terms of the number of places) groups of transformation units is carried out.The search is combinatorial, comparable to the MSS-method (cf.Sect.3.3), but with the difference that (1) the combined elements are the transformation units, (2) the consistency of the combined units is tested by a variance-analytical method, (3) only neighbouring units are combined.
The test of consistency of a sample of n of N given transformation units U k corresponds to a hypothesis test for a reduction of the number of parameters in the adjustment model (e.g.Niemeier, 2002, pp. 171).For the sample is tested, whether the 2n shifts Λ 0k and Φ 0k (k = 1...n) can be reduced to 2 shifts Λ 0 and Φ 0 .
From the n single adjustments of the U k a joint estimation s 2 0T of the variance of unit weight can be obtained: where s 2 0k and the redundance r k refer to U k .(Condition is that the s 2 0k are comparable, which is given here.)A second estimation s 2 0H = S H /r H of the variance of unit weight is obtained from an adjustment of the n regarded U k , in which joint shifts Λ 0 and Φ 0 are used (implicit hypothesis, see Niemeier, 2002, pp. 174).Of interest is the increase R = S H − S T of the weighted sum of squared residuals.Test statistic is T F will be F h,r T -distributed, if the involved U k have the same shifts.In practice, the computation of R is performed by a kind of an adjustment in steps (Niemeier, 1979, pp. 55), where the available results of the adjustments of the U k are used.
Condition for a merging of a sample of n transformation units U k is that the U k compose a continuous zone.To test this, methods of the graph theory are applied.An undirected graph is generated whose vertices stand for the U k .Between two vertices for U k and U l an edge will be set, if U k and U l lie geometrically neighboured.Then, the considered U k will be regarded as geometrically linked continuously, if the vertices for the U k form a connected graph.This is tested by means of the depth-first search (see, e.g.Turau, 2004, pp. 107).
For two transformation units U k and U l the condition of neighbourhood is tested by two tests: 1.For U k and U l polygons are generated which represent the units and their neighbourhood (cf.Sect.3.4).They are tested for overlap by (a) a test for intersection of their sides and (b) the point-in-polygon-test (cf.Sect.3.4), which covers the case that one polygon contains the other.
2. It is tested, whether the spherical distance between the centers of U k and U l is lower than a given maximal value d.
In practice, no mergings larger than n = 3 U k have occurred.Thus, to save computing time, the number of U k combined can be limited to, e.g.n = 5.

Postprocessing (step 6)
The in steps 2, 3, and 5 automatically generated transformation units are visualized and their formation is checked for improvements in case it is indicated.Above all, single places causing strong overlaps of transformation units or having large distances from the other places of their group are questionable results.The causative place is assigned to another transformation unit, and the new configuration is tested by a new adjustment and tests Eqs. ( 6), (7), and (8).If this is unsuccessful, also a scribal error can be quested on the basis of the estimated error Eq. ( 9).

Software
The analysis method described was implemented as a software package.The programming language is C. It is distinguished by the speed of its compiled programming code, which is useful for the computationally intensive combinatorially working algorithms used here.The results of the analysis can be visualized by various outputs: transformation units by means of DXF, KML, and command and data files for the Generic Mapping Tools; relative shifts of transformation units and residual vectors by means of DXF; transformed positions by means of KML.

Simulated example
The analysis method was tested by a simulated data set.The aim was to reconstruct given groups of places of homogeneous distortions and to detect wrong identifications and wrong ancient coordinates.Further, the accuracy of the transformation was determined.A data file with the data set can be found in the Supplement of this article.
www.hist-geo-space-sci.net/3/99/2012/ Hist.Geo Space Sci., 3, 99-112, 2012 The modern positions were transformed by Eq. ( 1) using the given parameters, and the resulting coordinates were contaminated by approximately normally distributed errors.The standard deviation of the errors in latitude was set to 5 and in longitude to 5 /cos45 • = 7 , where 45 • is the average latitude of the region considered.To degrade the conditions, the contaminated coordinate values were not rounded according to the highest precision occurring in the Geography but in longitude one third of the coordinates to the next 1 6 • , 1 4 • , and 1 3 • , respectively, and in latitude two thirds to the next 1 6 • and one third to the next 1 4 • .The resulting coordinates are the simulated ancient positions of the given places.In the following, a place is referred to by its number (no consecutive numeration; the hundred specifies the transformation unit).
The given transformation units, parameter values, and modern positions are regarded as the desired solution.
The generated ancient coordinates are regarded as original variants in the following.For 20 places (24 %) additional variants were generated (e.g. 5 instead of 30 , 0 instead of 20 ). 10 of these 20 places have > 1 identifications.The ancient coordinate variants can be regarded, for example, as values adulterated by scribal errors.Additionally, 4 an- The accuracy of the transformation of ancient into modern coordinates depends on the accuracy of the determined parameters.After the application of the analysis method to  the example, the accuracy of this transformation was determined.In doing so, ancient coordinates of additional places were generated randomly in the area of the given places.
Since only the accuracy of the transformation mattered in this test, the generated ancient coordinates were not contaminated.For a generated position a transformation unit was chosen by a point-in-polygon test based on the bufferpolygons of the given transformation units, where buffer width was b = 1 2 • (cf.Sect.3.4).(In the case of more than one possible units one was chosen randomly.)1000 pairs of positions and transformation units were generated.Then transformation Eq. (3) was performed using the desired parameters as well as the determined parameters.From the differences between both results the empirical variances were derived, where λ and φ are based on the desired parameters and λ and φ are based on the determined parameters.

Applying the analysis method
In the following, the results of the several runs of the analysis method are documented.Table 2, with a ≥ 1 6 • for Φ it is set σ Φ = 6.5 .That also holds for Λ; another approach, however, is σ Λ = σ Φ /cosφ ≈ 9 so that the intermediary value σ Λ = 8 was chosen.An adjustment of average parameters based on all given locations including the modern and ancient variants was performed.The result is: The deviations of the adjusted scales from the desired values are < 0.1.To degrade the conditions, for the following analysis the somewhat worse approximative values m (0) λ = 1.3 and m (0) φ = 1.0 were used.Further, to degrade the conditions, for the residual plot based on the adjustment result the additional and not the original ancient coordinate variants and about 50 % of the wrong modern positions were used.The residual plot and the constructed initial subsets of places are shown in Fig. 7 (modern position with ancient correction vector (v Λ v Φ ) ).Among the initial subsets there are four incorrect mergings of transformation units: E1 and E2, E4 and E5, E8 and E9, E10 and E11.

Modified MSS-method:
Variants of ancient coordinates cannot be considered by the method.To degrade the conditions, in the cases of ancient coordinate variants not the original but the additional variant was used.The parameters of the statistical tests are w max = 2.5 (lower value here for more reliable results), As a result, the desired transformation units formed partly, see Fig. 8; initial subsets containing incorrectly merged transformation units split up.20 unassigned locations remained.There were 1 wrong identification (No. 1349), 1 location in wrong transformation unit (No. 307 in E5 instead of E3), 4 accepted additional/wrong ancient coordinate values (Nos. 229,339,651,823).3. Forward-strategy: 19 locations were merged with transformation units, thereof No. 1313 with E12 instead of E13 due to its additional ancient coordinate variant.

Verification of the scales:
The adjusted parameters are m λ = 1.18 ± 0.02 and m φ = 1.09 ± 0.02.They do not deviate significantly from the initial values.

Merging of transformation units:
The parameter of the statistical test is α = 5 %; the parameters of the geometric tests are b = 1.5 • .According to Table 2, the new standard deviations are 14 and 9 , respectively.

Detection of scribal errors
The unassigned place No. 1314 was transformed (modern into ancient) by the parameters of E13, which came into consideration due to the position.The resulting ancient latitude is Φ = 48 • 37 ≈ 48 • 40 .In comparison with the given Φ = 47 • 40 the scribal error 48 • → 47 • became evident.Due to its larger σ Λ introduced in run 4, No. 1037 had a low influence on the adjustment so that the not detected scribal error in Λ is unproblematic.However, the resulting v Λ = 28 of Λ = 54 • 00 is noticeably large.The estimation of the error by Eq. ( 9) is ∇ = −28 so that the corrected Λ is 54 • 28 ≈ 54 • 30 and the scribal error 54 • 30 → 54 • 00 can be assumed.
The further two locations (Nos.823, 1312) with simulated scribal errors have no noticeable residuals.The reason is that the scribal errors correct the randomly adulterated coordinates to some extent.

Results
The nature of the residual plot of the example (Fig. 7) corresponds to that of the residual plots of the original data in Figs. 2, 3, 4.That argues for the kind of distortions assumed for the original data.
Despite the involved additional ancient coordinates and wrong identifications the scales were determined sufficiently accurately by a first adjustment of average values for the scales and shifts (step 1).The analysis method iterated within three runs to scale parameters being statistically identical with the given values.The automatic steps (2-5) of the method worked satisfactorily; except for one wrongly assigned place, the given groups of locations with homogenous shifts and the given modern identifications were reconstructed finally.By a manual preprocessing the wrong assignment was found and corrected.
There are three places with ancient coordinate variants whose given original variant was not determined by the analysis.However, the three determined variants are, compared to the original variants, more or comparably accurate (e.g.No. 339: given λ = 22 • 47 , systematically adulterated 47 • 20 , original (random) variant Λ = 47 • 30 , additional variant 47 • 10 ).Similarly, two of four simulated scribal errors were not detected, because the altered values are actually corrections.Thus, in practice, possibly not the original but the more accurate ancient coordinate variant is determined, which is preferable for a determination of accurate transformation parameters.
Table 3 gives the determined shifts in the form of relative shifts ∆ Λ0 and ∆ Φ0 with respect to the shifts Λ0E1 and Φ0E1 of E1, where ∆ Λ0 = Λ0 − Λ0E1 , ∆ Φ0 = Φ0 − Φ0E1 .The determined relative shifts agree well with the given values (Cols.∆Λ 0 , ∆Φ 0 ), the differences range from 0 to 5 .Furthermore, Table 3 gives the estimated standard deviations s Λ and s Φ of the ancient coordinates.In the case of transformation units with a larger number of places the accuracies s Λ = 6 to 12 and s Φ = 5 to 9 express the adulterations generated by coincidence and rounding well.In the case of less places the accuracy is overestimated.
The transformation of the 1000 randomly generated ancient locations by means of the determined transformation parameters resulted in accuracies of s λ = 2.3 and s φ = 2.1 by Eq. ( 16).Thus, the errors of the transformation are smaller than those of the ancient coordinates (here ≥ 5 ) so that they are negligible.

Conclusions
Since the geographic coordinates of the ancient places in Ptolemy's Geography are strongly erroneous, they cannot be used directly to find the modern counterparts of unknown Hist.Geo Space Sci., 3,[99][100][101][102][103][104][105][106][107][108][109][110][111][112]2012 www.hist-geo-space-sci.net/3/99/2012/ C. Marx: Rectification of the ancient geographic coordinates in Ptolemy's Geographike Hyphegesis 111 places.An insight into the distortions of the coordinates was given by a visualization for exemplary regions.For a rectification it is important to model the occurring errors close to reality.This is taken into account by a new analysis method which was presented.The method provides a way to validate the Ptolemaic coordinates of certainly identified places as well as to verify uncertain identifications.Further, by means of the determined systematic errors, ancient coordinates of so far unknown places can be transformed approximately into modern coordinates.In conjunction with topographic, historical, and archaeological information these calculational results can lead to the modern counterparts of the ancient places.
The analysis method was tested on a simulated example comparable with the real data situation.The method worked satisfactorily, the simulated systematic and crucial gross errors were determined.The transformation of ancient into modern coordinates turned out to be satisfactorily accurate.
In practise, the accuracy of the transformed coordinates is additionally downgraded by the (local) inaccuracy of the ancient coordinates considered.Further, deviations from the used error model or occasionally less accurate transformation parameters can adulterate the transformation result.Hence, the transformed coordinates and the derived identifications are afflicted with uncertainties.That is taken into account by giving a grade of certainty or the kind of information underlying the modern identifications.
The analysis method has been applied in an interdisciplinary project at the Technische Universität Berlin.The results for the identifications of the places of Book II are published in Kleineberg et al. (2010Kleineberg et al. ( , 2012)), the results for Book III will be published soon.
It is conceivable that the basic approach of the analysis method presented can be used or modified respectively for other problems concerning the partition of data into consistent subsets and gross erroneous elements, not only in view of historical data, but also today's measurement data.

2 3•
leads to σ = 2 3 • : 3 ≈ 14 .Standard deviations of other resolutions were interpolated by a logarithmic function, see σ(a) in Table2.a is not known in each case (

Figure 6 .
Figure 6.Simulated example: transformation units and relative shifts with respect to E1.

•Figure 7 .
Figure 7. Simulated example: residual plot and initial subsets of places.

Table 2 .
Model for the relation between resolution a and standard deviation σ of the ancient coordinates; σ(a) = d + log b (a + c), b = 1.214, c = 4.277, d = −7.508.