Matchmaker, Matchmaker, Make Me a Match: Migration of Populations via Marriages in the Past

The study of human mobility is both of fundamental importance and of great potential value. For example, it can be leveraged to facilitate efficient city planning and improve prevention strategies when faced with epidemics. The newfound wealth of rich sources of data—including banknote flows, mobile phone records, and transportation data—has led to an explosion of attempts to characterize modern human mobility. Unfortunately, the dearth of comparable historical data makes it much more difficult to study human mobility patterns from the past. In this paper, we present an analysis of long-term human migration, which is important for processes such as urbanization and the spread of ideas. We demonstrate that the data record from Korean family books (called “jokbo”) can be used to estimate migration patterns via marriages from the past 750 years. We apply two generative models of long-term human mobility to quantify the relevance of geographical information to human marriage records in the data, and we find that the wide variety in the geographical distributions of the clans poses interesting challenges for the direct application of these models. Using the different geographical distributions of clans, we quantify the “ergodicity” of clans in terms of how widely and uniformly they have spread across Korea, and we compare these results to those obtained using surname data from the Czech Republic. To examine population flow in more detail, we also construct and examine a population-flow network between regions. Based on the correlation between ergodicity and migration in Korea, we identify two different types of migration patterns: diffusive and convective. We expect the analysis of diffusive versus convective effects in population flows to be widely applicable to the study of mobility and migration patterns across different cultures.

Popular Summary

Physics-based approaches can be used to examine human migration. For instance, the movement of people over short and medium time scales can be traced in the recent past via phone records and the diffusion of bank notes. However, the historical mobility of humans over long time scales has remained nearly unexplored. We use quantitative tools based on physics to study historical human mobility and migration using an unusual data source: records from Korean family books (“jokbo”) that date back over 750 years. We develop notions for measuring the ergodicity of clans (i.e., the degree to which geographical diffusion has taken place) and are able to successfully identify and provide quantitative estimates for the extent of two distinct migration patterns.

Jokbo contain marriage records: When a marriage occurs, the bride’s clan name is recorded in the jokbo that belongs to the groom’s family. Because the extended version of Korean family names (e.g., “Kim from Gimhae” and “Lee from Hakseong,” to use examples of the names of researchers involved in this study) indicates a clan’s geographical origin, the data contained in jokbo can be used to trace not only the temporal history of clans but also their geographical movements. In our study, we use ten jokbo, which contain a total of over 220 000 entries that date back to the 13th century, as well as modern Korean census results from 1985 and 2000. We compute a “marriage flux”—the rate of marriages of women in some clan $i$ to men in some clan $j$—and we use these data to analyze population flows among 199 standardized regions of Korea. We find that some clans have reached a sort of “equilibrium” situation by being spread out geographically throughout Korea (3% of clans have at least one member in each region of Korea). The fraction of such ergodic clans increases with distance between the clan’s center of origin and Seoul, the crowded capital of Korea. To explain this result, we posit that the clans exhibit a bulk “convective” flow towards big cities and that clans have more opportunities to spread out geographically if they have to travel a further distance to Seoul. We also find that natural “diffusive” flow occurs. Our results suggest, however, that it would take longer than 50 000 years to produce a well-mixed Korean society by diffusion alone, so convective effects must also be at play because Korean culture has only existed for around 4000 years.

Our study incorporates both historical data, which tend to be both sparse and noisy, and modern (and well-sampled) census results. The use of historical data is important for investigations of long-term migration, which can inform studies of urbanization and the spread of ideas.

Figure 1

Examples of (a) ergodic and (b) nonergodic clans. We color the regions of South Korea based on the fraction of the total population composed of members of the clan in the year 2000. We use arrows to indicate the origins of the two clans: Gimhae on the left and Ulsan (“Hakseong” is the old name of the city) on the right. In this map, we use the 2010 administrative boundaries [39]. See the appendices for discussions of data sets and data cleaning.

Figure 2

Flux predictions from the population-product model (i.e., the special case of the gravity model with $\gamma =0$) with $\alpha =1$ and the radiation models for jokbo 1. (a) Scatter plot of the number of clan entries in jokbo 1 versus the corresponding centroid in 2000 using the population-product-model flux with $\alpha =1$. We compute the line using a linear regression to find the fitting parameter $a_G\approx 6.55(4)\times 10^{-11}$ (with a 95% confidence interval) to satisfy the expression $N_i=a_G{G}_{ij}$, where $G_{ij}$ is the population-product-model flux of women from clan $i$ to clan $j$ and $N_i$ is the total number of entries from clan $i$ in the jokbo. (b) We compare the same clan entries using the radiation model. We compute the line using a linear regression to find the fitting parameter $a_R\approx 0.049(2)$ to satisfy the expression $N_i=a_R{R}_{ij}$, where $R_{ij}$ is the radiation-model flux of women from clan $i$ to clan $j$ and $N_i$ is the total number of entries from clan $i$ in the jokbo. In both panels, we color the points using the number of administrative regions that are occupied by the corresponding clans [see Figs. 3 and 3]. The red markers (outliers) in both panels correspond to the clan of jokbo 1 (i.e., the case $i=j$).

Figure 3

Distribution of the number of different administrative regions occupied by clans. (a) Probability distribution of the number of different administrative regions occupied by a Korean clan in the year 2000. (b) Probability distribution of the number of different administrative regions occupied by the clan of a Korean individual selected uniformly at random in the year 2000. The difference between this panel and the previous one arises from the fact that clans with larger populations tend to occupy more administrative regions. [That is, we select a clan uniformly at random in panel (a), but we select an individual uniformly at random in panel (b).] Note that the rightmost bar has a height of 0.17, but we have truncated it for visual presentation. (c) Probability distribution of radii of gyration (in km) for clans in 2000. (d) Probability distribution of radii of gyration (in km) for clans of a Korean individual selected uniformly at random in 2000. The difference between this panel and the previous one arises from the fact that clans with larger populations tend to occupy more administrative regions. Solid curves are kernel density estimates (from Matlab R2011a’s ksdensity function with a Gaussian smoothing kernel of width 5).

Figure 4

Distribution of estimated diffusion constants (in ${\mathrm{km}}^2/\mathrm{year}$) computed using 1985 and 2000 census data and Eq. (3). The solid curve is a kernel density estimate (from Matlab R2011a’s ksdensity function with default smoothing). See Appendix pp7 for details of the calculation of diffusion constants.

Figure 5

Fraction of ergodic clans and distance scales of clans. For this figure, we use the 3 900 clans from the 2000 census data for which we were able to identify the origin location (see Appendix pp4). (a) Fraction of ergodic clans versus distance to Seoul. The correlation between the variables is positive and statistically significant. (The Pearson correlation coefficient is $r\approx 0.83$, and the $p$-value is $p\approx 0.0017$.) For the purpose of this calculation, we call a clan “ergodic” if it is present in at least 150 administrative regions. We estimate this fraction separately in each of 11 equally sized bins for the displayed range of distances. The gray regions give 95% confidence intervals. (b) Fraction of ergodic clans versus the distance between the location of clan origin and the present-day centroid. We measure ergodicity as in the left panel, and we estimate the fraction separately for each range of binned distances. (We use the same bins as in the left panel.) The correlation between the variables is positive and significant up to 150 km ($r\approx 0.94$, $p\approx 0.0098$) and is negative and significant for larger distances ($r\approx -0.98$, $p\approx 2.4\times 10^{-4}$).

Figure 6

Scatter plots of the number of clan entries in jokbo 2–10 versus the corresponding centroid in 2000 using the population-product-model flux with $\alpha =1$. We show our results in numerical order of the jokbo in panels (a)–(i), so jokbo 2 is in panel (a), etc. In each panel, we calculate the line using linear regression to determine the fitting parameter $a_G$ for $N_i=a_G{G}_{ij}$, where $G_{ij}$ is the population-product-model flux of women from clan $i$ to clan $j$ and $N_i$ is the total number of entries from clan $i$ in the given jokbo. The parameter values are (a) $a_G\approx 2.36(1)\times 10^{-9}$ [jokbo 2], (b) $a_G\approx 6.6(1)\times 10^{-11}$ [jokbo 3], (c) $a_G\approx 5.15(5)\times 10^{-9}$ [jokbo 4], (d) $a_G\approx 5.15(5)\times 10^{-8}$ [jokbo 5], (e) $a_G\approx 5.8(1)\times 10^{-9}$ [jokbo 6], (f) $a_G\approx 5.1(2)\times 10^{-16}$ [jokbo 7], (g) $a_G\approx 4.25(5)\times 10^{-8}$ [jokbo 8], (h) $a_G\approx 1.44(1)\times 10^{-9}$ [jokbo 9], and (i) $a_G\approx 4.71(8)\times 10^{-8}$ [jokbo 10]. The red markers in panels (a), (c), and (h) correspond to the clans of the depicted jokbo, and $N_i{|}_{i=j=\text{own}\text{clan}}=0$ for all of the other jokbo. In each case, we use a 95% confidence interval and color the points according to the number of administrative regions occupied by the corresponding clans.

Figure 7

Scatter plots of the number of clan entries in jokbo 2–10 versus the corresponding centroid in 2000 using the radiation-model flux. We show our results in numerical order of the jokbo in panels (a)–(i), so jokbo 2 is in panel (a), etc. In each panel, we calculate the line using a linear regression to determine the fitting parameter $a_R$ for $N_i=a_R{R}_{ij}$, where $R_{ij}$ is the radiation-model flux of women from clan $i$ to clan $j$ and $N_i$ is the total number of entries from clan $i$ in the jokbo. The parameter values are (a) $a_R\approx 0.062(2)$ [jokbo 2], (b) $a_R\approx 0.0098(7)$ [jokbo 3], (c) $a_R\approx 0.040(3)$ [jokbo 4], (d) $a_R\approx 0.075(7)$ [jokbo 5], (e) $a_R\approx 0.23(2)$ [jokbo 6], (f) $a_R\approx 0.0069(5)$ [jokbo 7], (g) $a_R\approx 0.12(1)$ [jokbo 8], (h) $a_R\approx 0.11(1)$ [jokbo 9], and (i) $a_R\approx 0.11(1)$ [jokbo 10]. The red markers in panels (a), (c), and (h) correspond to the clans of the depicted jokbo, and $N_i{|}_{i=j=\text{own}\text{clan}}=0$ for all of the other jokbo. In each case, we use a 95% confidence interval and color the points according to the number of administrative regions occupied by the corresponding clans.

Figure 8

For each jokbo, we plot the number of distinct clans $N_c$ versus the total number of entries $N_e$ on a doubly logarithmic scale. We calculate the red line via a linear regression to Heaps’ law [65] using the expression $N_c=10^{b_J}{N}_e^{a_J}$. This yields a slope of $a_J\approx 0.55(7)$ and an intercept of $b_J\approx 0.6(3)$ (with 95% confidence intervals).

Figure 9

Fraction of entries in each jokbo as a function of the birth year of brides using (a) linear and (b) semilogarithmic scales. The sudden drop on the right of each panel simply reflects the fact that women who are too young are not yet married.

Figure 10

Probability distribution for the number of different administrative regions occupied by the 822 “new” clans that are in the 2000 census data but are not in the 1985 census data. The solid curve is a kernel density estimate (from Matlab R2011a’s ksdensity function with a Gaussian smoothing kernel of width 1).

Figure 11

Probability distribution for how far clans have moved in terms of the distance from the historical clan origin location to the clan centroid from 2000. We geographically identified the origin and centroid for 3 900 clans among the 4 303 clans in the 2000 census data. The rightmost bar corresponds to all distances of at least 500 km, and the solid curve is a kernel density estimate (from Matlab R2011a’s ksdensity function with default smoothing).

Figure 12

Locations of clan origin names that we determined using the Google Maps API on top of the Korean map. We show administrative boundaries of South Korea that correspond to the upper-level local autonomy (광역자치단체) composed of provinces (도), a special autonomous province (특별자치도), a special city (특별시), and metropolitan cities (광역시) [73].

Figure 13

Distribution of clan-centroid movement between 1985 and 2000, which we compute using $\parallel {\mathbf{r}}_{i,\mathrm{C}}^{\mathrm{pop}}(t=1985)-{\mathbf{r}}_{i,\mathrm{C}}^{\mathrm{pop}}(t=2000)\parallel$ from Eq. (h1). (The norm $\parallel ·\parallel$ is the Euclidean norm.) The curve marked by squares weighs each clan equally, and the curve marked by circles weighs each clan by its population. We indicate distance in units of populations in 2000. (The 1985 data give a very similar distribution.) As one can see in the inset (for which we use a doubly logarithmic scale), the maximum movement distance is larger than 400 km. However, as the main panel illustrates, most clans moved considerably shorter distances.

Figure 14

We show box plots for the distributions of (a) $U_{\mathrm{tot},k\to k^{\prime }}$ for total population, (b) $U_{i,k\to k^{\prime }}$ for Kim from Gimhae, and (c) $U_{i,k\to k^{\prime }}$ for Lee from Hakseong, where the boxes indicate quartile values and the whiskers give the minimum and maximum values. We also show the probability distributions for (d) $U_{\mathrm{tot},k\to k^{\prime }}$ for total population, (e) $U_{i,k\to k^{\prime }}$ for Kim from Gimhae, and (f) $U_{i,k\to k^{\prime }}$ for Lee from Hakseong. In all cases, we apply cspline smoothing in gnuplot to the histograms using 50 equally-sized bins.

Figure 15

(a) Maximum relatedness subnetwork (MRS) [84] of the combined population-flow network [with adjacency-matrix elements (h6)] for all clans, (b) a magnified portion of this MRS that includes all regions with nonzero in-degrees, and (c) an instance of a model (inspired by the one in Ref. [64]) of a rewired version of a monocentric network (with Seoul as the center) with a rewiring probability of $p=0.4$. (See the discussion in the text.) Red lines indicate directed edges towards Seoul, and blue lines indicate directed edges towards other regions.

Figure 16

(a) MRS [84] of the population-flow network for (a) Kim from Gimhae and (b) Lee from Hakseong. In panel (a), red lines indicate directed edges towards Seoul, and blue lines indicate directed edges towards other regions. In panel (b), we show directed edges towards Ulsan in purple and directed edges towards other regions in blue.

Figure 17

Distributions of clan-density anomalies $\{{\varphi }_i\}$ in 2000 over the 199 standardized administrative regions for (a) Kim from Gimhae and (b) Lee from Hakseong. The leftmost and rightmost peaks correspond, respectively, to all values less than or equal to $-500$ and all values greater than or equal to 500. Solid curves are kernel density estimates (from Matlab R2011a’s ksdensity function with a Gaussian smoothing kernel of width 20).

Figure 18

Scatter plots of (a) distance between clan origin location and population centroid versus number of administrative regions, (b) distance between clan origin location and population centroid versus radius of gyration, and (c) radius of gyration versus number of administrative regions. (In this figure, we consider the 3 481 clans for which data are available in both the 1985 and 2000 censuses; we were able to determine the origin location for 3 120 of these clans.) The corresponding Pearson correlation values are (a) $r\approx -0.0094$ (from the 3 120 clans for which we know the origin location; the $p$-value is $p\approx 0.60$), (b) $r\approx 0.14$ (from the 3 120 clans for which we know the origin location; the $p$-value is $p\approx 1.1\times 10^{-15}$), and (c) $r\approx -0.26$ (from the 3 481 clans; the $p$-value is $p\approx 1.8\times 10^{-53}$). Note that correlations over limited ranges can be different and significant. For example, in panel (c), the two diagnostics for ergodicity are significantly positively correlated when $r_g<50\mathrm{km}$ ($r\approx 0.39$, $p\approx 2.3\times 10^{-4}$).

Figure 19

(a) Probability distribution of the number of different administrative regions occupied by a Czech family name in 2009. Note that the two leftmost bars have heights of 0.17 and 0.03 (with 0.17 at the farthest left), but we have truncated them for visual presentation. (These data were initially analyzed in Ref. [25].) (b) Probability distribution of the number of different administrative regions occupied by the clan of a Czech individual selected uniformly at random in 2009. The difference between this panel and the previous one arises from the fact that clans with larger populations tend to occupy more administrative regions. [That is, we select a clan uniformly at random in panel (a), but we select an individual uniformly at random in panel (b).] (c) Probability distribution of radii of gyration (in km) of Czech family names in 2009. Note that the leftmost bar has a height of 0.11, but we have truncated it for visual presentation. (d) Probability distribution of radii of gyration (in km) of Czech family names of a Czech individual selected uniformly at random in 2009. The difference between this panel and the previous one arises from the fact that clans with larger populations tend to occupy more administrative regions. Observe that the distributions in panels (a) and (b) are starkly different from the distributions in panels (a) and (b) from Fig. 3. Solid curves are kernel density estimates (from Matlab R2011a’s ksdensity function with a Gaussian smoothing kernel of width 5).

Figure 20

Scatter plots of diffusion constants that we compute using the 1985 and 2000 census data versus (a) distance between origin and population centroid, (b) number of administrative regions, and (c) radius of gyration of clans in 2000. (In this figure, we consider the 3 481 clans for which data are available in both the 1985 and 2000 censuses; we were able to determine the origin location for 3 120 of these clans.) The corresponding Pearson correlation values are (a) $r\approx -0.02$ (from the 3 120 clans for which we know the origin location; the $p$-value is $p\approx 0.21$), (b) $r\hspace{0.17em}\approx \hspace{0.17em}0.088$ (from the 3 481 clans; the $p$-value is $p\approx 1.7\times 10^{-7}$), and (c) $r\approx 0.097$ (from the 3 120 clans; the $p$-value is $p\approx 1.0\times 10^{-8}$).

Figure 21

Density-equalizing population cartograms [63] for South Korea using population data from (a) 1970 and (b) 2010 censuses [38]. The coordinates are longitude on the horizontal axis and latitude on the vertical axis. The growth of the Seoul metropolitan area over the past 40 years is clearly visible. (Compare this figure to a regular map of South Korea, such as the one in Fig. 1 in the main text.)

Figure 22

Radii of gyration and distance scales of clans. (For this figure, we use the 3 120 clans that are present in both the 1985 and 2000 censuses and for which we could determine the origin location.) (a) Radius of gyration $r_g$ versus distance to Seoul. The Pearson correlation between the variables is not statistically significant ($r\approx 0.18$; the $p$-value is $p\approx 0.6$). (b) Radius of gyration $r_g$ versus distance between the clan origin location and the present-day centroid. The Pearson correlation between the diagnostics is positive and statistically significant up to 170 km ($r\approx 0.86$, $p\approx 0.01$) and is negative and significant for larger distances ($r\approx -0.96$, $p\approx 0.005$). For each of the panels, we estimate $r_g$ separately in each of 11 equally-sized bins for the displayed range of distances. The gray regions give 95% confidence intervals.

Figure 23

Fraction of ergodic clans and distance scales of clans using only the 3 120 clans that we employed for the calculations in Fig. 22. We obtain the same qualitative result as in Fig. 5 in the main text. (a) Fraction of ergodic clans versus distance to Seoul. The correlation between the variables is positive and statistically significant. (The Pearson correlation coefficient is $r\approx 0.70$, and the $p$-value is $p\approx 0.02$.) For the purpose of this calculation, we call a clan “ergodic” if it is present in at least 150 administrative regions. (b) Fraction of ergodic clans versus the distance between the location of clan origin and the present-day centroid. We measure ergodicity as in the left panel, and we estimate the fraction separately for each range of binned distances. (We use the same bins as in the left panel.) The correlation between the variables is positive and significant up to 170 km ($r\approx 0.99$, $p\approx 0.0001$) and is negative and significant for larger distances ($r\approx -0.92$, $p\approx 0.01$). For each of the panels, we estimate the fraction of ergodic clans in each of 11 equally-sized bins for the displayed range of distances. The gray regions give 95% confidence intervals.