CSIS Discussion Paper No. 33
Analysis of Surface Changes by Tessellations
Yukio Sadahiro
NOVEMBER, 2000

Center for Spatial Information Science
University of Tokyo
7-3-1, Bunkyo-ku, Tokyo 113-8656, Japan

November 21, 2000
Analysis of Surface Changes by Tessellations

Yukio Sadahiro
Center for Spatial Information Science and Department of Urban Engineering
University of Tokyo
7-3-1, Hongo, Bunkyo-ku, Tokyo 113-8656, Japan
Phone: +81-3-5841-6273
Fax: +81-3-5841-8521
E-mail: sada@okabe.t.u-tokyo.ac.jp

1

Analysis of Surface Changes by Tessellations
Abstract
This paper develops a method for analyzing the change of surface, a continuousfunction defined over the two-dimensional space. The surface change is represented by ascalar field in a spatio-temporal region. The region is divided into subregions called hills
and waves which are based on the local form of the surface. They are both threedimensional tessellations of the spatio-temporal region, and provide an overview of thesurface change. To summarize the tessellations, four operations are proposed: section,
projection, fusion, and graph representation. They extract useful information fromtessellations to describe the structure of surface change. To test the validity of the method,
the change of a retail cluster in Shinjuku and Shibuya area in Tokyo is analyzed. The
empirical study yields some interesting findings that help us understand changes in thespatial structure of retailing.

2

1. INTRODUCTION
The change of spatial phenomena has drawn much attention of geographers.
Residential areas gradually expand with urban growth; retail stores gather around theaters
and hotels to form clusters and become shopping malls; people move from small to large
cities, which often leads to overconcentration to a few megacities; road networks growfrom simple, linear forms to more complicated, web-like forms. These changes have beenfrequently studied in geography and other related fields.
Analysis of changes, however, was rather impeded by lack of spatio-temporal data
in digital format. They were available only in limited cases, and consequently,
quantitative analysis in digital environment required manual digitization of time-seriesmaps. The cost of digitization was a major obstacle in spatio-temporal analysis by GIS.
Fortunately, this problem is now being solved by the progress of spatial data acquisition
tools such as GPS and mobile GIS (Longley et al. 1999). The cost of data acquisition is
far less expensive today and no longer an impediment to spatio-temporal analysis.
Improvement of data availability can promote spatio-temporal analysis in digital
environment. Its methodology, however, is still in the process of development. Most of
existing studies perform either spatial or temporal analysis fixing time or location, which
is not spatio-temporal analysis in its strict sense (Langran 1992). A summary statistic
representing spatial phenomena is calculated at discrete times or time periods and itschange is discussed on the temporal axis. Otherwise, a quantitative measure of temporalchange is calculated at sample points and its spatial distribution is examined. Thoughthese methods are useful in geographical analysis, there are some cases where theycannot describe the change of spatial phenomena. For example, basic transformations
such as translation and rotation are not detectable if we fix time in analysis. This isbecause spatio-temporality is not explicitly considered. It is obviously desirable to treatchanges in their original form, fixing neither space nor time.
This paper proposes a new method for analyzing the change of spatial phenomena,taking its spatio-temporality into consideration. Among various phenomena this paper
focuses on the surface, a scalar function defined over a two-dimensional region. Surfaces
are used for modelling population distribution (Tobler 1964; Griffith 1981; Bracken
1993; Bracken and Martin 1995), elevation of the earth's surface (Pike 1988; Hutchinson
1989; Etzelmuller 2000), distribution of geological measures (Isaaks and Srivastava
1989; Cressie 1993; Bailey and Gatrell 1995), and so forth. The surface is one of themost important objects in GIS, and its change is of interest to geographers.
To discuss the change of spatial phenomena, we have to consider a spatial region
and a time period simultaneously, that is, a spatio-temporal region. In Section 2, we
propose two types of tessellations of the spatio-temporal region where a surface function3

is defined. They are called hills and waves whose definition is based on the local form of
the surface. Hills and waves, however, are not easy to analyze directly, because they areboth three-dimensional objects. We thus propose four operations to summarize thetessellations in Section 3: section, projection, fusion, and graph representation. These
operations help us understand the structure of tessellations, and consequently, the changeof the surface. To test the validity of the method, Section 4 analyzes the change of a retailcluster in Shinjuku and Shibuya area in Tokyo, Japan. Section 5 summarizes theconclusions with discussion.
2. TESSELLATIONS OF SPATIO-TEMPORAL REGION
The surface is a continuous function defined over a two-dimensional region.
Consequently, if it changes over time, it is observed in a three-dimensional region, that
is, a spatio-temporal region. This leads to an analysis of the surface in a spatio-temporal
region.
To this end, we will propose two types of tessellations of the spatio-temporal
region where a surface function is defined. They are used to describe the surface change.
Prior to the proposal, however, we first introduce a tessellation of a two-dimensional
region, in order to help understanding of the tessellations discussed later.
Suppose a continuous function f(x) defined over a connected region S, which issmooth and twice differentiable at any point in S (Figure 1a). Approximating the local
surface around every point in S by a quadratic function, we obtain critical points where
df ( x )= 0. Critical points are classified into three categories: peaks, bottoms, and colsdx
(Warntz 1966; Okabe and Masuda 1984). From every point in S an integral curve of thesteepest ascent vectors of f(x) runs up to one of the peaks of f(x) (Figure 1b; Warntz
1966; Warntz and Waters 1975; Pfaltz 1976; Okabe and Masuda 1984). This naturally
gives a tessellation of S in which each subregion is assigned to one of the peaks of f(x)(Figure 1c).
FIG. 1. A surface and a tessellation of the region. (a) A surface function, (b) peaks (whitecircles) and integral curves of steepest ascent vectors, (c) a tessellation of the region.
Keeping this tessellation in mind, we proceed to tessellations of the spatio-temporal
region. Let S and T be regions defined in the spatial and temporal dimension,respectively. We consider their Cartesian product S×T and a continuous scalar function
f(x, t) in S×T which is twice differentiable at any point by both x and t. Differentiating
f(x, t) by x, we obtain a set of critical points discussed above. Since f(x, t) is a4

continuous function, they form a finite number of loci in S×T by which we can trace themovement of critical points. From all the loci we focus on those of peaks, and calculate
the tessellation of S×T as follows. We first fix t at t0, and assign every point in S to oneof the peaks to which an integral curve of the steepest ascent vectors of f(x, t0) runs up.
Moving t0 over T, we obtain an assignment of all the points in S×T. We finally reassign
every point in S×T to one of the loci to which its previously assigned peak belongs. Twopoints in S×T are assigned to the same locus if they were assigned to peaks in the same
locus. This procedure yields a tessellation of the spatio-temporal region S×T which wecall hills. In each hill all the points belong to the same locus. Figure 2 shows an exampleof hills where S is a one-dimensional region.
FIG. 2. Hills (regions bounded by dotted and solid lines) defined by f(x, t) in S×T ,where S is a one-dimensional region.
As shown in Figure 2, hills are bounded by two types of faces. A Boundary is aconnected set of points that divide S at every t which is depicted as the dotted line in
Figures 1 and 2. An Edge is another type of face that appears with generation anddisappearance of hills (solid lines in Figure 2). When a hill disappears, a set of points
around the peak are suddenly reassigned to another hill. This yields a face perpendicular
to the temporal axis.
Hills represent the history of two-dimensional subregions shown in Figure 1c.
Generation and disappearance of a subregion is indicated by edges of hills. Boundaries
tell us how the tessellation changes over time. Section of hills at a time gives theirsnapshot as a two-dimensional tessellation.
Hills are three-dimensional objects in the spatio-temporal region, though they areillustrated as two-dimensional in Figure 2. Similarly, both boundaries and edges are twodimensional objects. Consequently, hills and their boundaries do not always have simple
shapes; a hill may have another hill inside itself.
Interchanging S and T in the above discussion, we obtain another tessellation of thespatio-temporal region S×T. Differentiating f(x, t) by t, we obtain a set of temporal peaks∂f∂2 f= 0 and< 0 . Temporal peaks, as well as spatial peaks, form aof f(x, t) where∂t∂t 2
finite number of connected faces in S×T. Fixing x at x0, we assign every point in T to oneof the temporal peaks by climbing up the steepest ascent vectors of f(x0, t), which yields
a tessellation of the one-dimensional region T (Figure 3). Therefore, reassignment ofpoints to the faces to which the peaks belong yields a tessellation of S×T by the faces, theset of temporal peaks. We call the tessellation waves (Figure 4).

5

FIG. 3. Tessellation of the one-dimensional region T. (a) The section of f(x, t) at x=x0,(b) tessellation of T at x0 by f(x0, t).
FIG. 4. Waves defined by f(x, t) in S×T, where S is a one-dimensional region.
Waves have also two types of faces as shown in Figure 4. A connected set ofpoints dividing T at every x is called a boundary (dotted lines in Figure 3 and 4); todistinguish it from that of hills, we may say a wave boundary. The other is called an edge
or a wave edge, which appears at the spatial margin of waves. The wave edge is a concept
corresponding to the hill edge.
A wave is a connected set of points in S×T that share the change pattern of f(x, t). A
set of locations S0 belong to the same wave for a limited time period T0 if f(x, t) in S0changes similarly during T0. Let us consider, for example, spatial distribution of rent in abusiness district. When business is going well, rent rises everywhere in the district. In
depression phase, however, rent gradually goes down. Though there is a slight time lag
among locations, the business district shares the change pattern of rent. Therefore, thedistrict belongs to the same wave during the cycle of business.
Hills and waves are dual concepts of tessellations in the spatio-temporal dimension.
Consequently, they have many similarities such as boundaries and waves. However, thedifference of dimension between space and time causes several differences between them.
For example, the generator of a hill is a one-dimensional object, that is, a locus of spatial
peaks. On the other hand, a wave is assigned to a two-dimensional object - a faceconsisting of temporal peaks. Because of this, only waves can be divided further into
small pieces - former and latter halves, which often have an important meaning as shownin Section 4. In addition, waves can be arranged by their order on the temporal axis;
waves compose a semiordered set.
3. SUMMARY REPRESENTATION OF TESSELLATIONS
Hills and waves proposed in the previous section are both three-dimensional
objects in the spatio-temporal region. Unfortunately, it is not easy at present to analyze
three-dimensional objects directly because of the following reasons: 1) it is difficult tovisualize them in a three-dimensional space; 2) compared to two-dimensional objects,
three-dimensional objects convey a great amount of information; 3) GIS is not yet able todeal with three-dimensional objects in vector format. To avoid these difficulties, we
propose several methods for summarizing tessellations which help us understand theirstructure.

6

3.1 Section
A simple but useful method to analyze three-dimensional objects is to make theirsections by planes and lines. Fixing time t, we can visualize the structure of tessellations
on a plane. The section by line is usually made along the temporal axis at a fixed location.
If we add identification numbers to hills and waves, we can trace the surface change, say,the movement of hills, on a set of sections. Sections are effective at an early stage ofanalysis, especially in exploratory spatio-temporal analysis.
3.2 Projection
A more sophisticated method of summarization is to make projections oftessellations on S and T. A projection of hill boundaries on S indicates the extent of their
movement. A projection of hill edges on T shows the frequency distribution of generationand disappearance of hills.
Let us consider, for example, the projection of hill boundaries B on S. Let 1B(x, t)
be the indicator function of boundaries, that is,1 if a boundary exists at ( x, t ).
1B ( x, t ) = otherwise0(1)

An ordinary projection of B on S is given by
1 if ∃t ∈ T , 1B ( x, t ) ≠ 0.
P0 ( x, B; S ) = otherwise0(2)This projection represents the shadow of boundaries by the light parallel to the temporalaxis. This definition, however, is not useful if hills are so small that there are manyboundaries in S×T; the function P0(x, B; S) is equal to one almost everywhere.
To avoid this problem, we slightly modify the definition of projection operation:
objects are smoothed by a distance-decay weighting function in projection. The projectionof hill boundaries is then mathematically defined by
P( x, B; S ) = ∫ ∫ 1B (y, t )w( x − y )dydt ,(3)t ∈T y ∈S

where w(r) is a weighting function. One example of the weighting function is a Gaussianfunction r 2
w(r ) = exp −   .
(4)
  R 
The projection of boundaries on T is given by
P(t, B; T ) = ∫ ∫ 1B ( x, u)w' ( t − u )dudx ,x ∈S u ∈T(5)where w'(r) is a weighting function. Projections of edges can also be made in a similarway.

7Projections on S represent spatial stability of hills and waves. If hills do not change
their shape and location in T, the boundaries give dark shade in their projection, andconsequently, P(x, B; S) shows a wild fluctuation. Dark shade of wave edges indicatesthat waves occur at the same location.
Projections on T, on the other hand, reflect the distribution of topological events generation and disappearance of hills and waves. Concentration of topological events
makes the projection partially dark. If topological events are uniformly distributed in T,
the projection shows a uniform light shade.
3.3 Fusion
The above two methods reduce the amount of information conveyed bytessellations. These methods, however, do not necessarily keep important and useful
information about the surface change; it is difficult to make appropriate sections oftessellations before understanding their structure. One method to assure that is to use
fusion operation, that is, to combine small hills and waves into larger ones. This is basedon the assumption that larger hills and waves are more important than smaller ones inanalysis of surface changes.
In fusion operation, evaluation of the size of hills and waves is critical. A simple
but reasonable measure is their volume which is defined as follows. Let Hi be the ith hillin S×T, whose indicator function is given by
1 if ( x, t ) ∈ Hi.
(6)
1Hi ( x, t ) = 
0 otherwise
The volume of Hi is defined byvol( Hi ) = ∫∫ 1 (x, t )dxdt .

t ∈T x ∈S

Hi(7)Hills whose volume is less than a certain value are regarded as insignificant and combined
into neighboring hills.
In some applications, the surface height is also important. To take it into account,
we can use the height of Hi defined byheight ( Hi ) = max f ( x, t ) − min f ( x, t ),(8)
its absolute intensity

int ( Hi ) = ∫x ∈H ix ∈H i∫f ( x, t )1Hi ( x, t )dxdt ,t ∈T x ∈S

and its relative intensity
int ( Hi ) = ∫∫t ∈T x ∈S{( )}f x, t − min f ( x, t ) 1Hi ( x, t )dxdt .
x ∈H i(9)(10)

In addition to the above, there are several measures useful for distinguishingsignificant hills from insignificant ones: the area of projection of hills on S, the length ofprojection of hills on T, and so forth. The choice of the measure depends on the objective

8

of analysis.
Hills and waves judged to be small and insignificant are combined into larger ones.
There is a choice of hills into which small hills are combined. One reasonable choice isthe largest hill among neighboring ones. This method is used to eliminate spurious
polygons in GIS (Star and Estes 1990). If the surface height is important, we combine asmall hill into the neighboring one that gives the smallest relative depth of boundaries(Figure 5). This approach is often taken in topological analysis of surfaces and temporalfunctions (Okabe and Masuda 1984; Sadahiro 2000; Okabe and Masuyama 2001).
FIG. 5. Merging a small hill into larger one. H1, H2, H3: hills, B2, B3: boundaries. (a)Relative depth of boundaries with respect to the hill H1, (b) H1 is combined into H2
because the boundary B2 is shallower than B3.
3.4 Graph representation
If the topological structure of tessellations is of interest, graph representation ofhills and waves works effectively in analysis. Hills and waves are represented by nodes,
while boundaries and edges by links. For example, waves shown in Figure 2 arerepresented by the graph in Figure 6.
FIG. 6. Graph representing the topological structure of waves shown in Figure 2. Waves
are represented as nodes, while boundaries and edges by links.
Graphs are useful for not only description of the global structure of tessellations but
also comparison of tessellations. Let us consider two spatio-temporal regions in both ofwhich surface functions are defined. We say that the surface changes are topologicallyequivalent if the graphs representing hills and waves are isomorphic between the regions.
We can also evaluate the difference between surface changes by a distance measuredefined for comparison of graphs.
3.5 Combination of several methods
The methods proposed above can be used not only separately but also together inanalysis. Suppose that we are interested in the spatial structure of significant hills. To
remove insignificant hills, we apply the fusion operation to combine small hills into larger
ones. We then describe the topology of hills by a graph, and may finally make aprojection of the graph on S. If we are interested in the stability of hill boundaries, weagain remove insignificant hills by fusion operation and make a projection of hillboundaries on S. Combination of several methods extends applicability of hills and9waves.
4. EMPIRICAL STUDYThis section empirically analyzes the change of a surface in order to test the validityof the method proposed in Sections 2 and 3. We investigate the change of retail activity inShinjuku and Shibuya area in Tokyo, Japan (Figure 7). The study region has the biggest
retail cluster in Japan, surrounded by urban parks, residential areas and business districts.
The Yamanote Line and Meiji Street run through the region, and there are five major
railway stations: Shinjuku, Yoyogi, Harajuku, Shibuya, and Ebisu.
FIG. 7. Shinjuku and Shibuya area in Tokyo, Japan.
Spatial data used in the analysis are based on the list of retail stores and restaurants
in the NTT telephone directory which is published monthly. We extracted their addresses
in September every year from 1990 to 99 in ASCII format, and converted them into point
data by geocoding. Figure 8 shows the number of retail stores and restaurants in the
study region from 1990 to 99.
FIG. 8. The number of retail stores and restaurants, 1990-99.
From the point data we calculated the density distribution of stores at time i(i=1990, 91, ..., 99) using  x − y 2
ijf ( x, ti ) = ∑ exp −
(11)
 ,

R  
j



where yij is the locational vector of store j at time i and R is the distance decay parameter.
We tried three values for R, 100m, 150m, and 200m, but show only the result of
R=150m due to the limitation of space. From the set of density distributions we estimatedthe surface change from 1990 to 99 by linear interpolation at lattice points of resolution
10m.
Figure 9 shows the density distributions of stores in 1990 and 99. There were two
big subclusters in the study region, one at the northeast of Shinjuku, and the other at thewest of Shibuya. In addition, there were smaller subclusters at the east and west of
Shinjuku, the east of Harajuku and the southwest of Ebisu. These subclusters existed inboth 1990 and 1999, which implies the global stability of retail activity in the region.
FIG. 9. Density distribution functions of retail stores and restaurants, 1990 and 99.

10

We then calculated the tessellations of the region, 1990-99. To eliminate
insignificant hills, we performed fusion operation to combine hills whose relative heightis less than 150 into the highest neighboring hills. As a result, we obtained 90 hills and83 waves.
Figure 10 shows the distribution of absolute intensity of hills and waves defined byequation (9). Intensity of hills decreases linearly with the rank, while that of waves hasan S-shaped figure. This indicates that the intensity of waves has larger variance than thatof hills.
FIG. 10. Hills and waves arranged in order of the intensity.
Figure 11 shows the sections of hills in 1990 and 99, including the graphs
representing the topological structure of hills. There were 22 hills in 1990 along
Yamanote Line; there were several hills that cannot be recognized in Figure 9. Comparing
the figures we notice that two hills disappeared and three were newly generated since1990. The topology of hills represented by links had also changed, especially aroundShinjuku and Ebisu stations.
FIG. 11. Sections of hills and their graph representation, 1990 and 99.
Figure 11 suggests that the spatial structure of retailing had slightly changed at the
local scale. To analyze the change in more detail, we extracted stable hills that had existed
throughout the period; other hills were combined into neighboring ones by fusion
operation. Average and standard deviation of lifetime of hills are 4.19 and 3.18 years,
respectively.
Figure 12 shows the topology of stable hills. The figure indicates, against our
expectation, that the surface had considerably changed its spatial structure. There are only
11 stable hills, half of the hills that existed in 1990. They are recognizable in Figure 9,that is, relatively large hills. From this we conclude that the spatial structure of retailing
was stable at the global scale while it had drastically changed at the local scale.
FIG. 12. Stable hills and their topology.
In the above discussion we have focused on the topological structure of hills and its
stability over time. We then move to the stability of hills and waves with respect to theabsolute location in the spatio-temporal region. To this end, boundaries and edges of hills

11and waves were projected on S and T. In calculation of the projections, hills and waves
were weighted by the sectional intensity whose definition for hill boundaries is given by
intS ( Hi , t ) = ∫

x ∈S

{( )

}

f x, t − min f ( x, t ) 1Hi ( x, t )dx .
x ∈H i

The projection of hill boundaries is thenPS ( x, B; S ) = ∫ ∑ intS ( Hi , t )∫ 1Bi (y, t )w( x − y )dydt ,
t ∈T

y ∈S

i

(12)

(13)

where 1Bi(x, t) is the indicator function of the ith boundary. Figures 13 and 14 show the
projections on S and T, respectively.
FIG. 13. The projections of hill boundaries and wave edges on S. (a) Hill boundaries, (b)wave edges.
FIG. 14. The projections of hill edges and wave boundaries on T. (a) Hill edges, (b) wve
boundaries.
In Figure 13a we notice that the projection of hill boundaries clearly reflects thesurface distributions in 1990 and 99 shown in Figures 9 and 11. This supports our earlier
observation - global stability of retail activity in the region. In the south of the study
region, however, the projection shows relatively small values. The spatial structure of
retailing is rather unstable in that region.
Let us compare the projections of hill boundaries and wave edges on S in Figure
13. If they have similar shapes, it means that changes occur within hill boundaries;
correlation between spatial patterns of the surface and its temporal change. If they are
different, surface changes involve several hills. Concerning Figures 13a and 13b, weshould say that they are partly similar. In the north of the study region, both boundaries
and edges run from the west to the north east of Shinjuku station; linear pattern appears to
the south of Yoyogi station in both figures. Changes of the surface are restricted by its
spatial structure - they usually occur within hill boundaries. In the south, on the other
hand, boundaries and edges have quite different, even reverse patterns; most boundaries
run from southwest to northeast, while edges run from southeast to northwest. There islittle spatio-temporal correlation in the surface change in that region.
Let us turn to Figure 14. The figures are quite similar; hill edges and wave
boundaries had both gradually increased from 1990, reached the peak in 1995, anddecreased thereafter. The peak of hill edges in 1995 indicates drastic generation anddisappearance of hills, that is, great changes in the spatial structure of retailing. The
projection of wave boundaries, on the other hand, should be discussed in relation to thenumber of retail stores and restaurants shown in Figure 8. Stores had continuously

12decreased from 1990 to 1995, and gradually increased until 1998. This caused
disappearance of old waves and generation of new waves in 1995. From theseobservations, we say that the spatial structure of retailing drastically changes when thenumber of stores reaches the bottom - when the retail activity is in depression.
Generation and disappearance of hills are significant changes in the spatial structureof retailing. One clue to understand these phenomena is, as shown above, that they occur
when the retail activity is in depression. To analyze it in more detail, we examine thelocation of hill generations in relation to the size of waves. From 1990 to 99, there are 24
generations of hills which are concentrated in the 10 largest waves out of 83. Since thegeneration seems likely to occur in large waves, we performed the binomial test to
examine whether the hill generation occurs more frequently in larger waves. The null
hypothesis of random distribution was rejected at the significance level 1%, which
supports our hypothesis. Hills are generated by large waves, not by small waves. Since
the same result was obtained for the disappearance of hills, we conclude that the spatial
structure of retailing is more likely to change by large, long-term trends. In other words,
small and short-term trends are not enough to change the spatial structure of retailing.
We finally test whether the generation and disappearance of hills occur more
frequently in the former half of waves rather than their latter half (recall Section 2). The
binomial test, however, did not reject the null hypothesis of complete randomness. The
distribution of generation and disappearance is not biased with respect to the temporalpeaks of waves.
5. CONCLUSION
In this paper we have developed a method for analyzing surface changes bytessellations. The change of a surface is represented by a scalar field in a spatio-temporal
region. The spatio-temporal region is divided into hills and waves whose definition isbased on the local form of the surface. They are both three-dimensional objects
representing the spatio-temporal structure of the surface change. Unfortunately, hills andwaves are not easy to analyze directly because of the difficulty in visualization andunderstanding of a great amount of information they convey. We thus proposed four
operations to summarize the tessellations: section, projection, fusion, and graphrepresentation. They extract useful information from the tessellations to describe thestructure of surface change. To test the validity of the method, we analyzed the change ofa retail cluster in Shinjuku and Shibuya area in Tokyo, Japan. The empirical study yieldedsome interesting findings that help us understand the change of the spatial structure of
retailing.
In geographical analysis, visual analysis plays an important role not only at an early13

stage of analysis but also in mathematical modelling. Hills and waves, however, aredifficult to visualize as three-dimensional objects though they give us a great amount of
information about the structure of surface change. Consequently, visualization of hills
and waves is inevitably accompanied by dimension reduction, some of which aredescribed in Section 3. The choice of the reduction method heavily depends on theobjective of analysis. Users should carefully choose an appropriate method taking theirobjective into account.
The empirical study shown in Section 4 is at an early stage as location analysis: wetreated all the stores and restaurants equally, neglecting their category, size, customers,
and so forth. This is because we are interested in the spatial aspect of retailing rather thanits attribute distribution and spatial correlation between them. However, the study can befurther extended to take more factors into account in analysis. For example, the density
distributions of stores and the tessellations of the region are both obtained from the point
data of stores and restaurants. Therefore, hills and waves inherit the attributes of their
original points, and their spatial structures can be analyzed in relation to their non-spatial
properties by multivariate statistical methods including spatial autocorrelation coefficients.
Another way to treat diversity in stores and restaurants is to classify them by their
attributes and calculate surfaces for individual categories. Comparison of tessellations
would reveal the difference and relationship among surface changes, typical patterns of
the surface change, and so forth.
We finally discuss some limitations of our method for further research. First,
surfaces are usually estimated from a set of data measured at sample points. To convert
the discrete data into continuous form, we perform spatial interpolation such as spline(Spath 1995), kriging (Isaak and Srivastava 1989; Cressie 1993), and kernel estimation(Silverman 1986; Scott 1992). In addition, as done in the empirical analysis, the surfacechange is estimated from a set of surfaces of different times by temporal interpolation.
This implies that hills and waves calculated from the real data depend on the interpolation
method. One practical solution to this problem is to use fusion operation to eliminate
small hills and waves which are relatively unstable against the change of interpolation
methods. We can greatly reduce arbitrariness in the results of analysis by focusing only
on significant hills and waves. Another method to avoid the difficulty is to acquire spatial
data at a high resolution in both spatial and temporal dimensions. This is now a feasible
option because we have various tools for spatial data acquisition as discussed in the
introductory section. However, the ultimate solution is to develop spatio-temporalinterpolation techniques suitable for analysis of surface changes (Sadahiro 2000). New
methods that explicitly take the behavior of surfaces into account should be developed infuture research. Second, hills and waves give tessellations of the spatio-temporal region

14

which are based the first derivative of the surface function. In some applications,however, it may be necessary to discuss more detailed topography of surfaces, say,
aspect of slopes, ridge and course lines, and other concepts calculated from higher orderderivatives of the surface. Improvement of the ability of representing surfaces is an
important subject for future research.

15

LITERATURE CITED
Bailey, T. C. and A. C. Gatrell (1995). Interactive Spatial Data Analysis. London:
Taylor & Francis.
Bracken, I. (1993). "An Extensive Surface Model Database for Population-RelatedInformation." Environment and Planning B 20, 13-27.
Bracken, I. and D. Martin (1995). "Linkage of the 1981 and 1991 UK Censuses Using
Surface Modelling Concepts." Environment and Planning A 27, 379-390.
Cressie, N. (1993). Statistics for Spatial Data. New York: John Wiley.
Etzelmuller, B. (2000). "On the Quantification of Surface Changes using Grid-based
Digital Elevation Models (DEMs)." Transactions in GIS 4, 129-143.
Griffith, D. A. (1981). "Modeling Urban Population Density in a Multicentered City."
Journal of Urban Economics 9, 298-310.
Hutchinson, N. F. (1989). "A New Procedure for Gridding Elevation and Stream Line
Data with Automatic Removal of Spurious Pits." Journal of Hydrology 106, 211232.
Isaaks, E. H. and R. M. Srivastava (1989). An Introduction to Applied Geostatistics.
New York: Oxford University Press.
Langran, G. (1992). Time in Geographic Information Systems. London: Taylor &Francis.
Longley, A., M F. Goodchild and D. J. Maguire, D. J. (1999). Geographic Information
Systems: Principles, Techniques, Applications and Management. New York: John
Wiley.
Okabe, A. and S. Masuda (1984). "Qualitative Analysis of Two-dimensional Urban
Population Distribution in Japan." Geographical Analysis 16, 301-312.
Okabe, A. and A. Masuyama (2001). "An Exploratory Method for Qualitative Trend
Curve Analysis - A Robust Method against Low Quality Data." International Journalof Geographical Information Science, to appear.
Pike, R. J. (1988). "The Geometric Signature: Quantifying Landslide-Terrain Typesfrom Digital Elevation Models." Mathematical Geology 20, 491-511.
Pfaltz, J. L. (1976). "Surface Networks." Geographical Analysis 8, 77-93.
Sadahiro, Y. (2000). "Event-Based Analysis of Changes in Surfaces." CSIS Discussion
Paper Series 31, Center for Spatial Information Science, University of Tokyo, Japan(http://www.csis.u-tokyo.ac.jp/dp/31.pdf).
Scott, D. W. (1992). Multivariate Density Estimation. New York: John Wiley & Sons.
Silverman, B. W. (1986). Density Estimation. London: Chapman & Hall.
Spath, H. (1995). Two Dimensional Spline Interpolation Algorithm. Wellesley, MA: A
K Peters.

16

Star, J. and J. Estes (1990). Geographic Information Systems: An Introduction.
Englewood Cliffs, NJ: Prentice-Hall.
Tobler, W. (1964). "A Polynomial Representation of the Michigan Population." Papersof the Michigan Academy of Science, Arts, and Letters 49, 445-452.
Warntz, W. (1966). "The Topology of a Socioeconomic Terrain and Spatial Flows."
Papers of the Regional Science Association 17, 47-61.
Warntz, W. and N. Waters (1975). "Network Representations of Critical Elements ofPressure Surface." Geographical Review 65, 476-492.

17

(a)

(b)

Figure 1(c)Tboundaryedge
locus of spatial peaks

S

Figure 2

T

T

(a)

(b)

Figure 3

Tboundaryedge
face consisting of temporal peaks

S

Figure 4

relative depth ofboundary B2

relativeve depth ofboundary B3
H2

H1
B2

H3

H2

H3

B3

B3

(a)

(b)

Figure 5Tboundary
edge
graph

S

Figure 6

Shinjuku

Yoyogi

Yamanote
Line

Meiji
StreetHarajuku
railway
Yamanote Linerailway
Meiji Street
study region
urban parks

Shibuya

Ebisu

Figure 754000

52000

50000

48000

46000

0
90

91

92

93

94

95

96

Figure 8

97

98

99

Shinjuku

Shinjuku

Yoyogi

Yoyogi

Harajuku

Harajuku

ShibuyaShibuya

Ebisu

Ebisu

0100 200 300 400 500 600 700 800 -

100
200
300
400
500
600
700
800

railway
road

1990

1999

Figure 9

1.0E+9

1.0E+91.0E+81.0E+81.0E+6

1.0E+6

1.0E+4

1.0E+4

1.0E+2

1.0E+2

1.0E+01.0E+0
1

20

40

60

80

90

1

20

40

60

rank

rank

intensity of hills

intensity of waves

Figure 10

80

90

peak
link
boundary
railway
road

1990

1999

Figure 11

Shinjuku
Yoyogi

Harajuku

Shibuya

Ebisu

Figure 12

peak
link
railwayroadShinjuku
Shinjuku
Yoyogi

Yoyogi

Harajuku

Harajuku

ShibuyaShibuya

Ebisu

Ebisu

(a)

(b)

Figure 13

1.6E+9

1.6E+91.2E+91.2E+90.8E+90.8E+90.4E+90.4E+9

0

0
90

91

92

93

94

95

96

97

98

99

90

(a)

91

92

93

94

(b)

Figure 14

95

96

97

98

99