Coastal Engineering 160 (2020) 103719

Contents lists available at ScienceDirect

Coastal Engineering
journal homepage: http://www.elsevier.com/locate/coastaleng

Probabilistic tsunami hazard assessment with simulation-based
response surfaces
T. Kotani a, K. Tozato b, S. Takase c, S. Moriguchi d, *, K. Terada d, Y. Fukutani e, Y. Otake b,
K. Nojima a, M. Sakuraba a, Y. Choe faResearch and Development Center, Nippon Koei Co., Ltd. Inarihara, 2304, Tsukuba-shi, Ibaraki, 300-1259, Japan
Department of Civil Engineering, Tohoku University, Aza-Aoba 468-1, Aramaki, Aoba-ku, Sendai, 980-0845, Japan
Department of Civil Engineering and Architecture, Hachinohe Institute of Technology, 88-1 Ohbiraki, Myo, Hachinohe, Aomori, 031-8501, JapandInternational Research Institute of Disaster Science, Tohoku University, Aza-Aoba 468-1, Aramaki, Aoba-ku, Sendai, 980-0845, JapaneCollege of Science and Engineering, Kanto Gakuin University, Mutsuura Higashi 1-50-1, Kanazawa-ku, Yokohama-shi, Kanagawa, 236-8501, Japanf
Department of Industrial and Systems Engineering, University of Washington, 3900 E Stevens Way NE, Seattle, WA, 98195, USAbc

A R T I C L E I N F O

A B S T R A C T

Keywords:
Probabilistic Tsunami hazard analysis
Response surface
Numerical simulations
Uncertainty analysis
Fault modelThis paper proposes a novel response surface-based method for probabilistic tsunami hazard assessment (PTHA).
Although recent advancements in numerical simulation have enabled the accurate characterization of tsunami
hazards, the high computational cost of numerical simulation often prohibits its broad application to probabilistic hazard analysis. The proposed method addresses this challenge by constructing the response surface (RS) of
a target output using the results of high-fidelity tsunami simulations. The proposed method quantifies uncertainties in key simulation variables and propagates those uncertainties to the target output through the RS in acomputationally efficient Monte Carlo simulation (MCS). We illustrate and validate the proposed method
through a case study of the tsunami induced by the 2011 Great East Japan earthquake. The case study focuses on
Sendai, Ishinomaki, and Kamaishi in Japan as the target locations. The proposed method estimates coastal
tsunami heights while considering uncertainties in the fault slip and rake as well as the modeling error associatedwith the numerical simulation. The MCS allows us to estimate the probability density functions of the tsunami
height at the target locations. The proposed method quantifies the contribution of each source of uncertainty tothe overall uncertainty in the target output and thus facilitates engineering decision-making.

1. Introduction
Preparing for catastrophic natural hazards is a significant engineering challenge. The hardest part of this challenge is striking a balance
between economy and safety; ensuring complete safety by structural
measures incurs a prohibitive cost. To address this challenge, probabilistic approaches to evaluating such hazards and risks have been
extensively discussed in the literature. In the field of seismic engineering, Cornell (1968) initiated an approach known as probabilistic seismic
hazard analysis (PSHA). Subsequently, this approach has been studied
by numerous researchers (e.g., McGuire, 1977; Anderson and Trifunac,
1978; Ishikawa and Kameda, 1988) and summarized in several reports
and books (e.g., National Research Council, 1988; Working Group on
California Earthquake Probabilities, 1995, 2002; McGuire, 2004).
Probabilistic tsunami hazard assessment (PTHA), which is based on

PSHA, enables us to quantitatively evaluate the probability of tsunami
risks with a certain return period; it is one of the practical tools being
used for disaster mitigation. Many studies on PTHA have been reported
in the literature (e.g., Power et al., 2007; Thio et al., 2007; Burbidge
et al., 2008; Mitsoudis et al., 2012; Suppasri et al., 2012; Stefanakis
et al., 2014; Fukutani et al., 2015; Park and Cox, 2016; Griffin et al.,2016; Park et al., 2018; Park et al., 2019) and summarized in some review papers (e.g., Geist and Parsons, 2006; Behrens and Dias, 2015;
Grezio et al., 2017; Mori et al., 2017). In particular, Annaka et al. (2007)
proposed a very simple practical method that has been implemented
successfully in PTHA handbooks currently used in Japan (Japan Societyof Civil Engineers (JSCE) 2011; JNES 2014). In this method, uncertainties are classified into two types: aleatory uncertainties due tovariations in actual phenomena and epistemic uncertainties due to insufficiencies in the knowledge, data, and modeling accuracy. Aleatory* Corresponding author.
E-mail address: s_mori@irides.tohoku.ac.jp (S. Moriguchi).
https://doi.org/10.1016/j.coastaleng.2020.103719
Received 25 June 2019; Received in revised form 30 March 2020; Accepted 28 April 2020
Available online 20 May 20200378-3839/© 2020 Elsevier B.V. All rights reserved.

T. Kotani et al.

Coastal Engineering 160 (2020) 103719

uncertainties are estimated based on the results of case studies andnumerical simulations, whereas epistemic uncertainties are often
handled using a logic-tree approach (e.g. Annaka et al., 2007; Burbidge
et al., 2008, Omira et al., 2016). The design standard published by theAmerican Society of Civil Engineers (ASCE) is particularly important inthis field. This standard (American Society of Civil Engineers, 2017)provides an extensive PTHA framework that is widely recognized as thefirst systematic PTHA framework. This framework (detailed in Sections
6.7 and C6.7 of the standard) provides a comprehensive guideline on
determining inundation depths and flow velocities based on site-specific
PTHA. This framework covers a broad range of uncertainties: epistemic
uncertainties in model parameters (e.g., magnitude, fault depth, fault
geometry, location, and slip distribution) and aleatory uncertainties (e.
g., natural variability in source processes, modeling uncertainties, andtidal variation). By contrast, for illustration purposes, this paper’s PTHA
focuses on the uncertainties in two fault parameters (i.e., slip and rake)
and the modeling error for the 2011 Great East Japan earthquake.
Geographically, the ASCE standard covers Alaska, Hawaii, and the Pacific Coast of the U.S. concerning tsunami loads and effects. The standard recommends the use of 1) a statistically weighted logic tree
approach to account for epistemic uncertainties and 2) probability distributions to account for aleatory uncertainties. For both types of uncertainties, this study uses probability distributions, as discussed later.
While the standard provides detailed guidance on numerical simulations
for PTHA (e.g., model spatial resolution), it does not directly advise on
computational costs of PTHA. This paper addresses this gap by proposing a surrogate model-based computationally efficient way to perform
PTHA.
Numerical simulations have also become essential tools for tsunami
hazard assessment. Currently, these simulations are highly accurate
thanks to recent advancements in computational mechanics (e.g., Larsen
and Fuhrman, 2019; Sarfaraz and Park, 2016; Qin et al., 2018). However, numerical simulations are currently not being fully utilized in
PTHA. Specifically, in the Japanese PTHA handbooks, the use of numerical simulations is restricted to obtaining the median value of theprobability density function of the maximum tsunami height. The meanand standard deviation of the probability density function are estimatedusing a different empirical method proposed in the literature (e.g., Aida,
1978). However, considering recent advancements in tsunami simulations, we postulate that numerical simulations can be further used forPTHA. Hence, we propose a new method that effectively leverages numerical simulations for PTHA.
The novelty of the proposed method for PTHA lies in its use of the
response surface (RS) of a target output, which is constructed based on
the results of highly accurate tsunami simulations. Because tsunami
simulations are computationally expensive, their direct use forcomprehensive PTHA is computationally prohibitive. The current
practice typically discretizes continuous probability distributions. For
example, simulation runs are performed at 100 randomly chosen
epicenter locations in a continuous geographic area (Mitsoudis et al.,2012); alternatively, probability masses are assigned to 72 tsunami
scenarios to sparsely represent a continuous spectrum of scenarios (Parkand Cox, 2016). This discretization makes PTHA computationallyfeasible at the expense of an approximation error. We address thischallenge by using the RS as a surrogate model of the tsunami simulation
model. The key idea is that the RS allows the continuous simulation
output to be accurately approximated by interpolating/extrapolating
outputs obtained at nearby inputs. Thus, once the important uncertainties are quantified (using empirical data, simulations, and engineering knowledge), an extensive Monte Carlo simulation (MCS) can be
performed with the RS to stochastically evaluate possible hazards andrisks. In our case study, we use only 50 numerical simulation runs toconstruct the RS, which is in turn evaluated 10,000 times with an MCS to
perform comprehensive PTHA.
To demonstrate its capability, the proposed method is applied to acase study of the tsunami induced by the 2011 Great East Japan

earthquake. Our PTHA considers uncertainties in the fault slip and rake
along with the modeling error in the numerical simulations. Sendai,
Ishinomaki, and Kamaishi in Japan are selected as target study locations.
The RS at each location is obtained based on numerical simulations ofcoastal tsunami heights. Finally, the probability density functions of thetsunami height are estimated by performing an extensive MCS.
Furthermore, we quantify the contribution of each source of uncertainty
to the overall uncertainty to determine the relative importance of eachsource for engineering decision-making. However, we note that the risk
assessment exemplified in this study may not be directly applicable toother practical cases where additional sources of uncertainties should betaken into account. A main objective of this paper is to demonstrate thespecific procedure of the proposed method and to discuss its capability
and potential for practical applications.
2. Methodology
2.1. Fundamental scheme of the proposed method
In general, comprehensively evaluating the random field of a target
output requires numerous simulation runs. Thus, using high-fidelity
numerical simulations for this purpose is computationally expensive, if
not infeasible. To address this problem, Honjo (2011) proposed ascheme, specifically for reliability-based design in geotechnical engineering. In this scheme, numerical simulations are carried out to obtaina RS that describes the trend of the target output subject to specificuncertainties. This scheme is highly efficient for performing a practicalreliability assessment and is the basis of the PTHA method proposed inthis study.
Fig. 1 shows the flowchart of the scheme proposed by Honjo (2011).
The scheme consists of three phases: (i) geotechnical analysis, (ii) uncertainty analysis, and (iii) reliability assessment. In the geotechnicalanalysis phase, a series of simulations are performed for different
calculation conditions and sets of input parameters. The results enable
the construction of the RS in a later phase; this RS models the relationship between the target output Y and the set X of key variables suchas important calculation conditions and input parameters. This process
informs the analyst of how each variable in X influences the target
output. The analyst also develops a sense of the physics and mechanisms
governing the output Y by analyzing the simulation results based on
his/her engineering knowledge and experience. In the uncertainty analysis phase, the uncertainties underlying X are quantified using databasesand engineering knowledge. For example, if X is a set of quantitative
variables, their probability distributions or general tendencies (e.g.,
average and standard deviation) should be determined. In contrast, if X

Fig. 1. Flowchart of the reliability-based design scheme proposed by
Honjo (2011).
2

T. Kotani et al.

Coastal Engineering 160 (2020) 103719

is a set of qualitative variables, the probabilities of observing qualitative
events are required. If the uncertainties in X are all captured in databases, then quantifying the uncertainties in X is straightforward; however, in practice, it is usually necessary to rely on engineering judgments
because of insufficient data. Finally, in the reliability assessment phase,
the RS of the target output Y is constructed using the results of thegeotechnical analysis. The RS is defined as a function of the key variables in X that describes the tendency of Y observed in the numericalsimulations. Using the RS with the quantified uncertainties from theprevious phase, a Monte Carlo simulation (MCS) is performed to obtain
the probability density function of the target output Y. Because the RS isa surrogate model of the output, which is constructed based on the
numerical simulations, the MCS results are expected to describe thesame tendency observed in the numerical simulation results. Note thatnumerical simulation runs should be designed in careful considerationof the range of X because the RS aims to model the tendency of Y within,
not beyond, the range of X.

tsunami simulations would be interesting, this error is not taken into
account in the present study. On the other hand, the modeling error
captures the discrepancy between an empirical phenomenon and its
numerical model; in other words, the modeling error represents theincompleteness of the model due to simplification and idealization andany inaccuracy due to discretization. This error is quantified in thisstudy by comparing the simulated results with observational data.
The uncertainty analysis can be performed in parallel with thetsunami simulations. The analysis begins by gathering statistical information of the key variables. This information (e.g., mean, variance, anddistribution type) is used to determine the probability distributions of
the key variables.
Subsequently, the reliability assessment is performed based on the
results of the tsunami simulations and uncertainty analysis. Using thetsunami simulation results, an RS is constructed as an approximate
function Y ’ FðXÞ. Note that the effects of bathymetry and topography
are represented directly in the tsunami simulations and hence are also
reflected in the RS. Then, we perform an MCS for the reliability
assessment using the RS and the probability distributions of the key
variables.
Some related studies in the literature (e.g. Sarri et al., 2012; Srajet al., 2014; Omira et al., 2015; Lorito et al., 2015; Selva et al., 2016;
Volpe et al., 2019) also statistically analyze simulation results to reduce
the total computational cost of PTHA. These studies, however, do notconsider the modeling error of numerical simulations. To the best of ourknowledge, this study is the first to quantify the effect of the modelingerror in comparison to other uncertainties.

2.2. Proposed PTHA method
In the reliability-based design scheme outlined in Section 2.1, thegeotechnical analysis phase is separate from the other two phases. In thisstudy, we propose replacing the geotechnical analysis phase with thetsunami analysis phase.
Fig. 2 shows the flowchart of the PTHA method proposed in thisstudy. The proposed method first specifies the target output. The output
can be quantitative such as the tsunami height, or qualitative, such aswhether an area is flooded. Then, the method specifies the key variables,such as important input parameters and calculation conditions used in
the tsunami analysis. The key variables should be determined based on
numerical simulations. Specifically, we can vary variables that are expected to significantly influence the target output based on engineeringknowledge and experience. It is efficient to determine the ranges of thesekey variables while observing how the variables affect the simulation
outputs. Note that a large number of simulation runs is not necessary
because the goal is simply to identify the tendency of the target output
with respect to key variables.
There are two important uncertainties in numerical simulations,
namely, transformation error and modeling error, according to Honjo(2011). The transformation error occurs when observational data aretransformed into the input parameters of numerical simulations because
the transformation is based on engineering assumptions and/or empirical knowledge. Although a discussion of the transformation error in

3. Application
3.1. Target output
In this section, a case study is presented to illustrate the specificprocedure of the proposed RS-based PTHA method and validate its
capability. This case study concerns the tsunami induced by the 2011
Great East Japan earthquake. We specified the tsunami height as thetarget output, and the proposed PTHA method was used to evaluate thetsunami heights at three coastal locations in Japan—Kamaishi, Sendai,
and Ishinomaki (see Fig. 3). These locations are close to the coastline
with water depths of 10–20 m. While this study focuses on the tsunami
heights along the coastline as the target outputs for illustration purposes,
future work can also use the proposed method to study tsunami runups.

Fig. 2. Flowchart of the proposed RS-based PTHA method.
3

T. Kotani et al.

Coastal Engineering 160 (2020) 103719

Fig. 3. Measurement points for the tsunami heights: (a) Kamaishi, (b) Sendai, and Ishinomaki.

3.2. Step (i): Tsunami analysis
3.2.1. Uncertainty in the tsunami hazard
Various phenomena characterize a tsunami (e.g., initial uplift of thesea surface induced by fault movement, followed by propagation andrunup), each with its own uncertainty. Because it is impractical toincorporate all possible uncertainties into a model, it is necessary to
select and model the key uncertainties. The fault model is recognized as
the most widely used model in tsunami analysis (Kajiura, 1963, 1982).
Various fault parameters, including the fault slip and rake, are illustrated in Fig. 4. A JSCE committee summarized the influence of eachfault parameter on the tsunami height (JapanSociety of Civil Engineers,
2002). Based on previous studies, we selected the fault slip and rake as
the two key variables in this study. To calculate the fault displacement,
we used the Fujii–Satake model Ver. 8.0 (FS55v80) (Satake et al., 2013),
which is one of the most widely used fault models for the 2011 Great
East Japan earthquake (see Fig. 5). This model considers 55 subfaults
and provides the spatial slip distribution. According to Goda et al.
(2014), this model’s tsunami predictions “agree well with the tsunami
observations.”

Fig. 5. Fujii–Satake model Ver. 8.0 (FS55v80) and the result of a centroid
moment tensor (CMT) inversion (after Satake et al., 2013).

3.2.2. Calculation conditions and governing equations
As explained in the previous subsection, the Fujii–Satake model Ver.
8.0 (FS55v80) considers 55 subfaults. We varied the slip and rake values
in the tsunami simulation to assess their effects on the target output. As
summarized in Table 1, we considered 50 test cases based on 5 slip
values and 10 rake values to cover reasonable ranges of the parameters.
The slip and rake values used in case S3R5 are the same as the original
parameters of FS55v80. The other calculation conditions are given in
Table 2. The initial distribution of the tsunami height relative to thestandard sea level in Japan (Tokyo Peil) was obtained using the method
proposed by Okada (1992), in which the vertical sea-surface displacement is calculated using a fault model. To describe tsunami propagation,
the following continuity equation (1) and nonlinear shallow waterequations (2) and (3) were used:∂η ∂M ∂Nþ
þ¼0
∂t ∂x ∂yFig. 4. Illustration of the fault parameters.
4(1)

T. Kotani et al.

Coastal Engineering 160 (2020) 103719

with grid intervals of 2430, 810, 270, and 90 m. For the area covered
with grids of intervals 2430–270 m, ETOPO1 digital elevation data
(Amante and Eakins, 2009) were used for the topographical conditions.
For the area covered with a grid interval of 90 m, we used digital
elevation data from the Digital 10 m Grid Elevation Map published by
the Geospatial Information Authority of Japan. We also used the M7000
series of isobathic data (JapanHydrographic Association, 2015), which
are high-resolution seafloor topography data. The reflecting boundary
condition was used on the coastline.

Table 1
Calculation cases..

3.2.3. Simulation results
We began by checking the validity of the tsunami analysis model
used in this study. As explained in Section 3.2.2, the conditions ofcalculation case S3R5 correspond to those of FS55v80. It is therefore
possible to compare the S3R5 results with the empirical data of the 2011
Great East Japan earthquake. The time histories of the wave heights at
the five offshore points indicated in Fig. 5 were used as the observationaldata, as shown in Fig. 7. These data were collected by GPS wave recorders of the Nationwide Ocean Wave Information Network for Portsand Harbours (The Port and Airport Research Institute, 2011). As seen in
Fig. 7, the simulation results generally match the observational data at
the five offshore points. Any discrepancies represent the modeling error
discussed in Section 2.2. Section 3.3.1 will detail how we quantified themodeling error. Note that because observational data are subject tomeasurement errors, the modeling error is also affected by measurement
errors.
The spatial distributions of the maximum tsunami height were
determined within 180 min after tsunami generation, and the results areshown in Fig. 8. Additionally, the maximum tsunami heights at Sendai,
Ishinomaki, and Kamaishi are summarized in Table 3 for the 50 cases.
These results indicate that fault slip is the dominant parameter at all
locations.

Table 2
Calculation conditions.
Item

Note

Scenario

Tsunami hazard of the 2011 Great East Japan
earthquake
Fujii–Satake model Ver. 8.0
Calculate crustal deformation and sea-surface
displacement using the method of Okada (1992)
Nonlinear shallow water equations
Leapfrog method and staggered grids
Domain 1: 2430 m; Domain 2: 810 m;
Domain 3: 270 m; Domain 4: 90 m
Offshore boundary: transmitting boundary
Landward boundary: reflecting boundary
Set at Tokyo Peil (T.P.) %0.57 m (coseismic tidal
level)
180 min from earthquake generation

Fault model
Initial waveform
Governing equations
Calculation scheme
Computational grid interval(with nesting)Boundary conditions
Tidal level
Calculation time

!

"

!

"

∂M ∂ M 2∂ MN∂η gn2 pﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃ
þ gD þ 7 M M 2 þ N 2 ¼ 0
þþ∂y D∂x D 3∂t ∂x D="!
! "∂N ∂ MN∂ N2∂η gn2 pﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃþþþ gD þ 7 N M 2 þ N 2 ¼ 0∂y D∂y D 3∂t ∂x D(2)

3.3. Step (ii): Uncertainty analysis
3.3.1. Quantification of uncertainties in key variables
Uncertainties in key variables, represented as probability distributions, are important inputs for the MCS. To specify the probability distributions, we can use statistical information such as the average μ,
standard deviation σ, and probability distribution type. There are two
approaches to obtain this information: 1) estimation based on empirical
data and 2) estimation based on engineering judgment and knowledge.
In this study, we employed the former approach for the modelingerror and the latter approach for the fault parameters (slip and rake).
Table 4 shows the estimated statistical values for the key variables.
In this study, normal distributions were used to represent epistemic(3)=

where η is the vertical displacement of the water surface; D is the total
water depth; g is the acceleration due to gravity; n is the Manning coefficient; and M and N are the fluxes in the x and y directions, respectively. The governing equations were discretized using the staggered
leapfrog scheme (Goto and Ogawa, 1982; Goto et al., 1997). The
calculation grids were generated as shown in Fig. 6 according to Imamura et al. (2006). The computational domain was divided into calculation grids that were nested from the wave source to the target locations

Fig. 6. Nested analysis regions.
5T. Kotani et al.

Coastal Engineering 160 (2020) 103719

Fig. 7. Comparison between the observational data and simulation results at the offshore points shown in Fig. 5 (Observational data: Nationwide Ocean Wave
Information Network for Ports and Harbors (The Port and Airport Research Institute, 2011).

uncertainties in the fault parameters (slip and rake). While the best
distribution for the slip may be right heavy-tailed according to an indepth study on the 2011 Great East Japan earthquake (Goda et al.,2014), our use of a normal distribution serves to illustrate how to
implement the proposed procedure when there is little to no empiricalknowledge of the fault parameters. The choice of a normal distribution isjustified according to the principle of maximum entropy when we canspecify only the mean and variance of a continuous distribution torepresent epistemic uncertainty (Dowson and Wragg, 1973).
As mentioned in Section 3.2.2, the slip and rake values of eachsubfault were expressed as ratios with respect to the original FS55V80
parameters. Thus, the values were normalized (see Table 1), and theaverage slip and rake values were set as 1.0. Based on the moment
magnitude scale Mw , the slip and rake standard deviations were specified as 0.1 and 0.04, respectively. According to Kanamori (1977), Mw isrelated to the slip D as follows:
log10 M0 ¼ 1:5Mw þ 9:1;

(4)

M0 ¼ rDS;(5)

used to cover this range.
The modeling error of the numerical simulation was also considered
one of the uncertainties. The time histories shown in Fig. 7 were used toquantify the modeling error Merr, which is defined as the difference
between the observed tsunami height Hobs and the simulated tsunami
height Hsim:
Merr ¼ Hobs % Hsim :(6)

Because we possessed time histories from both simulations and observations, the modeling error was quantified at numerous times. Fig. 9shows a histogram of the modeling error showing that the modelingerror can be represented as a normal distribution with an average andstandard deviation of approximately 0.49 and 0.76 m, respectively. Notethat there are some limitations on the modeling error considered in thisstudy; we considered only the time histories of the tsunami height at
offshore points and did not take other effects into account. For instance,
the rigid boundary conditions along the coastline have effects thatshould be considered. Other available information, such as tsunami
runups and inundation areas, may be considered in future work. Additionally, note that in this study, a normal distribution was assumed toroughly capture the tendency of the modeling error. However, in general, when choosing distributions to describe uncertainties, users of the
proposed procedure would want to consider various parametric distributions (e.g., a skew-normal distribution to better describe the positive
skewness observed in the modeling error depicted in Fig. 9) as well as
nonparametric distributions if the data size is sufficiently large.

where M0 is the seismic moment, r is the rigidity, and S is the fault area.
Hanks and Kanamori (1979) proposed a modified relationship, but theoriginal version, which is still widely accepted, is used in this study. As
reported by JapanSociety of Civil Engineers (2011), the estimatedstandard deviation of Mw of an earthquake having the same magnitude
as the 2011 Great East Japan earthquake is approximately 0.1. Using thisvalue, the range of values of the normalized slip becomes 0.71–1.41.
Because the range of a normal distribution is approximately &3σ , a slip
standard deviation of 0.1 was used to cover the range of 0.71–1.41. The
rake standard deviation was determined similarly. JapanSociety of Civil
Engineers (2011) conducted calculations with different rake values andreported that &10' is a suitable rake range; this range becomes 0.9–1.1
upon normalization. Therefore, a rake standard deviation of 0.04 was

3.4. Step (iii): Reliability assessment
3.4.1. Construction of the response surface
The RS of the tsunami height was constructed using the tsunami
simulation results in Table 3. In this study, the slip and rake wereselected as the independent variables of the RS. The RSs at the three6T. Kotani et al.

Coastal Engineering 160 (2020) 103719

Fig. 8. Spatial distributions of the maximum tsunami height for each calculation: (a) 90-m mesh shown in Fig. 6(a) and (b) 90-m mesh shown in Fig. 6(b).

target locations (Sendai, Ishinomaki, and Kamaishi) were constructed
separately because their responses were sufficiently different from eachother. Prior to constructing the RSs, we investigated the tendency ofeach response with respect to each independent variable. From thisinvestigation, we confirmed a linear relationship between the tsunami
height and the slip and a nonlinear relationship between the tsunami
height and the rake. These tendencies agree with the findings in relevant
studies (Necmioglu and Ozel, 2014; Goda et al., 2014). Therefore, we
added a quadratic term in the rake to the RSs. Additionally, a slip–rakecross term was included. The resulting form of the RS is
Hmax ¼ aU þ bλ þ cUλ þ dλ2 þ e;

RS efficiently for numerous independent variables. However, because of
the relatively modest number of terms considered in this study, we
checked every functional form using multiregression analysis.
Tables 5–7 show the fitted value of each coefficient, the adjusted R2
value, the residual standard error (r.s.e.), and the value of the Akaike
information criterion (AIC) (Akaike, 1973) of the RS for each of the threetarget locations in Fig. 3. We selected the functional form with thelowest AIC value as the most suitable one. As a result, the same functional form was selected coincidentally for Sendai, Ishinomaki, andKamaishi. The selected forms are highlighted in boldface in Tables 5–7.
The fact that the same functional form was selected for all the target
locations means that the dominant terms were identical independent oflocation. The constructed RS for each target location is shown in Fig. 10,
where the blue surface is constructed from the red points obtained fromthe numerical simulations (see Table 3). Because a visual inspection ofFig. 10 confirms that the variation in each RS is small between adjacent
red points, the number of simulation runs seems sufficient. Note that
determining the optimal design of computational experiments to build asurrogate model such as an RS is an active research topic in statistics (e.
g., Stefanakis et al., 2014).

(7)

where Hmax is the tsunami height, U is the slip, λ is the rake, and a, b, c, d,and e are undetermined coefficients. Note that the modeling error is notincluded in Eq. (7). The coefficients were determined by fitting the results in Table 3. There are four terms that include the independent RS
variables, and we considered 15 forms of RSs (all 24 forms except theconstant RS).
In general, we can use forward selection, backward elimination,
stepwise methods, or other variable selection methods to optimize the7

T. Kotani et al.

Coastal Engineering 160 (2020) 103719

Table 3
Maximum tsunami heights for the 50 cases at Sendai, Ishinomaki, and Kamaishi.
Case

Slip

Rake

S1R1
S1R2
S1R3
S1R4
S1R5
S1R6
S1R7
S1R8
S1R9
S2R10
S2R1
S2R2
S2R3
S2R4
S2R5
S2R6
S2R7
S2R8
S2R9
S2R10
S3R1
S3R2
S3R3
S3R4
S3R5
S3R6
S3R7
S3R8
S3R9
S3R10
S4R1
S4R2
S4R3
S4R4
S4R5
S4R6
S4R7
S4R8
S4R9
S4R10
S5R1
S5R2
S5R3
S5R4
S5R5
S5R6
S5R7
S5R8
S5R9
S5R10

0.700
0.700
0.700
0.700
0.700
0.700
0.700
0.700
0.700
0.700
0.850
0.850
0.850
0.850
0.850
0.850
0.850
0.850
0.850
0.850
1.000
1.000
1.000
1.000
1.000
1.000
1.000
1.000
1.000
1.000
1.200
1.200
1.200
1.200
1.200
1.200
1.200
1.200
1.200
1.200
1.400
1.400
1.400
1.400
1.400
1.400
1.400
1.400
1.400
1.400

0.735
0.815
0.877
0.938
1.000
1.062
1.124
1.185
1.247
1.309
0.735
0.815
0.877
0.938
1.000
1.062
1.124
1.185
1.247
1.309
0.735
0.815
0.877
0.938
1.000
1.062
1.124
1.185
1.247
1.309
0.735
0.815
0.877
0.938
1.000
1.062
1.124
1.185
1.247
1.309
0.735
0.815
0.877
0.938
1.000
1.062
1.124
1.185
1.247
1.309

Tsunami height [m]SendaiIshinomaki

Kamaishi

7.695
7.880
7.911
7.808
7.759
7.590
7.496
7.266
6.974
6.745
9.296
9.309
9.255
9.273
9.272
8.986
8.702
8.481
8.285
7.877
10.58710.80510.65810.61210.65410.281
10.205
9.859
9.370
9.107
12.40212.73512.732
12.66212.60312.331
11.926
11.51311.11910.51114.45414.41014.581
14.58214.43414.28713.87513.325
12.57512.041

4.226
4.443
4.626
4.780
4.897
4.982
5.033
5.049
5.034
4.981
5.036
5.295
5.516
5.694
5.834
5.927
5.982
5.994
5.975
5.919
5.824
6.125
6.384
6.598
6.762
6.878
6.943
6.954
6.907
6.824
6.827
7.185
7.487
7.744
7.942
8.085
8.171
8.211
8.169
8.064
7.830
8.226
8.552
8.823
9.057
9.226
9.338
9.384
9.366
9.263

10.003
9.950
9.845
9.693
9.484
9.210
8.951
8.593
8.197
7.714
11.80311.724
11.571
11.34511.07710.75710.369
9.962
9.537
9.046
13.47613.382
13.21012.95612.68112.30811.828
11.33410.737
10.16015.77015.672
15.46815.16114.836
14.39613.82613.257
12.537
11.77518.28218.09917.82717.51617.10916.52015.90615.20414.43813.553Fig. 9. Histogram of the modeling error approximated by a normal
distribution.

The modeling error is considered in each RS by adding Eq. (6) to Eq.
(7). Hence, the final form of the RS is as follows:
Hmax ¼ aU þ bλ þ cUλ þ dλ2 þ e þ Merr

The resulting RSs are used for the MCS in the next section.

3.4.2. Monte Carlo simulation
The MCS was performed using the RSs and the probability distributions of the key variables from the previous sections. Specifically, we
first generated the values of the key variables (slip, rake, and modelingerror) from the probability distributions. Using these values as inputs to
the RSs, we evaluated the tsunami height at each target location. Using
10,000 replications, we then estimated the probability density functionof the tsunami height.
3.4.3. Probability density of the tsunami heightFig. 11 shows the estimated probability density function of thetsunami height at each target location. Using the density function, we
can estimate the exceedance probability for any tsunami height. In thisstudy, as a validation of the proposed method, the tsunami heightobserved at each target location for the 2011 Great East Japan earthquake was used to estimate the exceedance probability; hence, thetheoretical exceedance probability should be 50%.
The tsunami heights observed at the target locations are summarizedin Table 8. These heights were derived from tsunami trace height datacollected in a field survey (Mori and Takahashi, 2012), from which weused the tsunami trace height reported at the survey point closest to eachtarget location except for Sendai; because three survey points werelocated near the target location for Sendai, the average tsunami traceheight was used. The tsunami heights estimated from the field survey
and the calculated exceedance probability are shown in Fig. 11. For
Sendai and Ishinomaki, the median values are close to the results of the
field survey. For Kamaishi, the median value seems much larger than thefield survey result (the exceedance probability is 91%). However, wenote that their absolute difference is 1.5 m, which is only approximately
10–15% of the field survey result. Because our numerical simulation
used the perfect reflecting boundary condition despite the complex
shape of the coastline near Kamaishi, we attribute this relatively large
deviation to a wave amplification effect near the coastline.

Table 4
Statistical information on the fault slip, rake, and modeling error.
Basic variables

AveragevalueStandard
deviation

Distribution function
type

Slip
Rake
Modeling
error

1.0
1.0
0.487

0.1
0.04
0.758

Normal distribution
Normal distribution
Normal distribution(8)

The results show that the slip affects the tsunami height more than
the rake does. This is in accordance with the findings reported in theliterature. Although the tsunami height decreases with increasing rake
at both Sendai and Kamaishi, the opposite tendency is seen at Ishinomaki. We reason that this difference is due to the positional relationshipbetween the locations of cities and the directions of tsunami propagation. Because the direction of the coastline at Ishinomaki is different
from those at the other two cities, the rake would have a different effect.

3.4.4. Calculation of contribution ratios
The variance of the probability density function of a target output
can be apportioned between multiple sources of uncertainties. It is
practically important to analyze the sensitivity of the target output withrespect to each uncertainty. In this section, the contribution ratio Ri ofeach uncertainty (i.e., in the slip, rake, and modeling error) is quantified
as follows:
8T. Kotani et al.

Coastal Engineering 160 (2020) 103719

Table 5
Response surface for Sendai.
No.

Coefficients of each term
a

b

C

d

e

20.2928 (14.665)
8.1858 (1.380)¡3.3048 (¡8.225)
8.4496 (27.649)%1.5426 (%1.692)¡9.6583 (¡15.124)%9.6583 (%3.371)%0.7472 (%1.600)%9.6583 (%9.664)¡9.1908 (¡11.561)
3.6591 (1.223)
1.9786 (3.653)%5.6817 (%5.412)
0.7693 (0.708)
7.7667 (26.03)
2.8170 (12.62)
4.2784 (12.45)
13.6193 (26.34)
3.659 (0.295)
1.1611 (6.365)
12.1578 (11.789)
1.1611 (3.001)
13.619 (7.032)
5.5320 (6.086)

R2

r.s.e. [m]

AIC

0.997
0.944
0.984
0.993
0.983
0.942
0.984
0.980
0.931
0.026
0.984
0.037
0.926
0.033
0.392

0.125
0.561
0.297
0.196
0.305
0.566
0.303
0.337
0.620
2.330
0.302
2.317
0.640
2.321
1.841

¡59.3089.92
26.41%15.42
28.98
89.95
27.43
38.09
98.96
231.3727.12
229.87101.35230.07206.87(t-value)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15

12.4757 (29.684)
10.6590 (11.162)
9.0688 (81.347)
12.4757 (12.169)

16.8888 (8.173)
0.3800 (0.366)

9.0688 (52.51)
9.0688 (47.23)%3.3048 (%3.372)
8.4719 (%21.45)%3.0240 (%11.25)%11.7271 (%19.40)
16.889 (0.686)

12.1162 (41.648)
9.0688 (24.837)%5.7217 (27.50)%1.5135 (%12.94)

8.4496 (25.03)%9.658 (%0.812)%2.9562 (%12.995)%3.024 (%1.633)%1.5135 (%1.693)

4.6807 (5.707)

Table 6
Response surface for Ishinomaki.
No.

Coefficients of each term
a

b

c

D

e

12.9900 (25.20)
9.14613 (4.811)

1.8460 (12.33)
5.57792 (56.970)
2.9741 (5.416)¡6.2871 (¡26.42)%6.28714 (%6.848)%0.5829 (%2.071)%6.28714 (%12.77)¡5.7844 (¡19.53)%1.70469 (%1.778)
1.3654 (4.184)%7.74449 (%14.97)
0.69916 (1.048)
2.8848 (25.11)%0.25095 (%1.503)%1.2609 (%6.08)
4.7789 (22.67)%1.705 (%0.213)
0.7276 (6.507)
5.7888 (8.657)
0.7276 (2.992)
4.779 (3.821)
2.1471 (7.04)

R2

r.s.e. [m]AIC

0.99910.98610.98620.996
0.98490.97950.97780.9822
0.97250.02497
0.98520.02743
0.93
0.03168
0.8356

0.04658
0.17970.1791
0.096410.18720.218
0.2267s0.20330.25271.504
0.18521.502
0.40291.499
0.6177¡158.03%23.90%24.27%86.19%19.81%5.52%1.60%12.51
9.24
187.62%21.81
186.5554.94
186.3397.68(t-value)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15

3.9609 (25.30)
2.7980 (4.864)
5.86388 (106.74)
3.96092 (6.293)

14.89138 (14.62)
0.02763 (0.043)

5.86388 (45.393)
5.8639 (50.62)

1.9290 (11.90)%3.8163 (%15.49)
14.891 (0.937)

3.9348 (22.062)
5.8639 (25.543)

1.929 (1.613)

1.84598 (3.068)
5.6029 (47.24)
5.5779 (40.52)
1.8713 (13.419)%1.8887 (%18.39)
0.89442 (10.227)%6.287 (%0.818)
0.8944 (1.543)

4.3514 (15.81)

Table 7
Response surface for Kamaishi.
No.

Coefficients of each term
a

b

C

d

e

18.4973 (21.21)
2.1523 (0.274)¡6.2884 (¡24.83)
9.5804 (23.676)%4.6821 (%5.907)¡8.8271 (¡21.93)%8.8271 (%2.326)%0.7045 (%1.735)%8.8271 (%5.783)¡7.5201 (¡15.01)
9.8278 (2.481)
2.6610 (5.650)%0.8431 (%0.526)
1.5829 (1.668)
10.9078 (28.15)
5.2057 (20.84)
8.2599 (20.84)
18.9307 (29.20)
9.828 (0.691)
1.8903 (11.86)
15.876 (13.481)
1.8903 (2.648)
18.931 (8.555)
7.865 (6.380)

R2

r.s.e. [m]AIC

0.99920.932
0.99180.989
0.99130.93330.98580.98140.92560.1226
0.99140.13750.82770.13320.23220.078810.74280.25840.29880.26650.73550.33910.38850.77692.669
0.2639
2.646
1.182
2.652
2.496

¡105.44344
117.99572
12.41913
26.94284
15.46987
116.07724
38.64625
52.26004
121.55742
244.95227
13.58715
243.14541
162.60447
243.39644
237.33548(t-value)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15

16.8426 (63.57)
15.1866 (18.293)
10.3600 (60.840)
16.8426 (18.803)
10.3600 (53.62)
10.3600 (46.80)
16.5605 (65.17)
10.3600 (15.377)

σ2( 100 ½%*;
σ all

Ri ¼ 2i

12.0202 (3.808)
0.2982 (0.329)%6.1789 (%19.94)%16.0467 (%21.18)
12.020 (0.427)%6.179 (%2.921)%6.2884 (%7.344)
9.5863 (23.96)
9.5804 (22.64)%6.0148 (%30.27)%7.7920 (%22.49)%3.0302 (%23.17)%8.827 (%0.648)%3.030 (%2.969)

4.423 (3.977)

on the uncertainty i and σ2all is the variance based on all the uncertainties.
The results obtained from Eq. (9) for each target location are given in
Table 9. For Sendai and Kamaishi, the contribution ratio of slip isapproximately 60%, and that of the modeling error is approximately(9)

where σ2i is the variance of the probability density function based only
9

T. Kotani et al.

Coastal Engineering 160 (2020) 103719

Fig. 10. Response surfaces.

Fig. 11. Exceedance probabilities.

30%. These results can be understood in light of the fact that, as
mentioned in Section 3.4.1, the fault slip is the dominant parameter forthe tsunami height. In contrast, for Ishinomaki, the contribution ratio ofslip is approximately 30%, while that of the modeling error is approximately 60%. This indicates that the inaccuracy of the numerical simulation at this target location contributes more to the overall uncertainty
than that at the other two target locations. A possible reason for this is
the effect of wave amplification near the coastline. However, to better
understand this phenomenon, further investigation is needed because
Kamaishi has an even more complex coastline shape than Ishinomaki.

The contribution ratios of the rake for all the target locations are
relatively small. This finding is also important because the uncertainty
in the rake can be less emphasized in engineering decision-making.
4. Conclusion
We proposed a novel method for the response surface-based probabilistic assessment of tsunami hazards that efficiently leverages
advanced numerical simulations. The method was illustrated and validated through a case study evaluating the tsunami induced by the 201110

T. Kotani et al.

Coastal Engineering 160 (2020) 103719

CRediT authorship contribution statement

Table 8Tsunami trace heights at the target locations (2011 Tohoku Earthquake Tsunami
Joint Survey Group, 2011).
Tsunami trace height [m]SendaiIshinomaki
Kamaishi

12.12911.2039.317
6.902
11.426

T. Kotani: Writing - original draft, Methodology. K. Tozato: Methodology. S. Takase: Methodology. S. Moriguchi: Conceptualization,
Writing - review & editing. K. Terada: Supervision, Writing - review &editing. Y. Fukutani: Methodology. Y. Otake: Methodology. K.
Nojima: Visualization. M. Sakuraba: Visualization. Y. Choe: Writing review & editing.

Position of tsunami trace
LatitudeLongitude

38.2183
38.2187
38.2184
38.413939.2650140.986140.986140.985141.287
141.892

Acknowledgments
The authors would like to gratefully acknowledge the support madethrough the University of Washington-Tohoku University: Academic
Open Space (UW-TU:AOS).

Table 9Contribution of each source of uncertainty to the target output uncertainty.
Targetpoint

Uncertainties (Keyvariables)

Varianceσ2[m2]

Contribution ratio[%]

Sendai

All
Slip
Rake
Modeling error
All
Slip
Rake
Modeling error
All
Slip
Rake
Modeling error

1.488
0.904
0.008
0.580
0.925
0.352
0.008
0.577
1.779
1.150
0.049
0.575

–
60.7
0.6
39.0
–
37.6
0.9
61.5
–
64.8
2.7
32.3

Ishinomaki

Kamaishi

References
Aida, I., 1978. Reliability of a tsunami source model derived from fault parameters.
J. Phys. Earth 26, 57–73. https://doi.org/10.4294/jpe1952.26.57.
Akaike, H., 1973. Information theory and an extension of the maximum likelihood
principle. In: Petrov, B.N., Caski, F. (Eds.), Proceedings of the 2nd International
Symposium on Information Theory. Akadimiai Kiado, Budapest, pp. 267–281.
Amante, C., Eakins, B.W., 2009. ETOPO11 arc-minute global relief model: procedures,
data sources and analysis. NOAA Technical Memorandum NESDIS NGDC-24. https://data.nodc.noaa.gov/cgi-bin/iso?id¼gov.noaa.ngdc.mgg.dem:316.
American Society of Civil Engineers, 2017. ASCE standard "Minimum design loads forbuildings and other structures" (ASCE/SEI 7-16). https://doi.org/10.1061/9780
784414248.
Annaka, T., Satake, K., Sakakiyama, T., Yanagisawa, K., Shutom, N., 2007. Logic-tree
approach for probabilistic tsunami hazard analysis and its applications to theJapanese coasts. In: Birkhauser Basel (Ed.), Tsunami and Its Hazards in the Indianand Pacific Oceans, pp. 577–592. https://doi.org/10.1007/s00024-006-0174-3.
Anderson, J.G., Trifunac, M.D., 1978. Uniform risk functionals for characterization ofstrong earthquake ground motion. Bull. Seismol. Soc. Am. 68 (1), 205–218.
Behrens, D., Dias, F., 2015. New computational methods in tsunami science. Phil. Trans.
R. Soc. A (Math. Phys. Eng. Sci.) 373. https://doi.org/10.1098/rsta.2014.0382,
2053.
Burbidge, D., Cummins, P.R., Mleczko, R., Thio, H.K., 2008. A probabilistic tsunami
hazard assessment for western Australia. Pure Appl. Geophys. 165, 2059–2088.
https://doi.org/10.1007/s00024-008-0421-x.
Cornell, C.A., 1968. Engineering seismic risk analysis. Bull. Seismol. Soc. Am. 58 (5),1583–1606.
Dowson, D., Wragg, A., 1973. Maximum-entropy distributions having prescribed firstand second moments (corresp.). IEEE Trans. Inf. Theor. 19 (5), 689–693.
Fukutani, Y., Suppasri, A., Imamura, F., 2015. In: Stochastic Analysis and Uncertainty
Assessment of Tsunami Wave Height Using a Random Source Parameter Model that
Targets a Tohoku-type Earthquake Fault, vol. 29. Stochastic Environmental Researchand Risk Assessment, pp. 1763–1779.
Geist, E.L., Parsons, T., 2006. Probabilistic analysis of tsunami hazards. Nat. Hazards 37,
277–314.
Goda, K., Mai, P.M., Yasuda, T., Mori, N., 2014. Sensitivity of tsunami wave profiles andinundation simulations to earthquake slip and fault geometry for the 2011 Tohoku
earthquake. Earth Planets Space 66 (105). https://doi.org/10.1186/1880-5981-66105.
Goto, C., Ogawa, Y., 1982. Tsunami Numerical Simulation with Leap-Frog Scheme.
Tohoku University, p. 52 (In Japanese).
Goto, C., Ogawa, Y., Shuto, N., Imamura, F., 1997. Numerical Method of Tsunami
Simulation with the Leap-Frog Scheme (IUGG/IOC Time Project), IOC Manual.
UNESCO.
Grezio, A., Babeyko, A.Y., Baptista, M.A., Behrens, J., Costa, A., Davies, G., Geist, E.L.,
Glimsdal, S., Gonz!
alez, F.I., Griffin, J., Harbitz, C.B., LeVeque, R.J., Lorito, S.,
Løvholt, F., Omira, R., Mueller, C., Paris, R., Parsons, T., Polet, J., Power, W.,
Selva, J., Sørensen, M.B., Thio, H.K., 2017. Probabilistic tsunami hazard analysis:
multiple sources and global applications. Rev. Geophys. 55 (4), 1158–1198. https://doi.org/10.1002/2017RG000579.
Griffin, J., Pranantyo, I.R., Kongko, W., Hauhan, A., Roniana, R., Miller, V., Davies, G.,
Horspool, N., Maemunah, I., Widjaja, W.B., Natawisjaja, D.H., Latief, H., 2016. In:
Assessing Tsunami Hazard Using Heterogeneous Slip Models in the Mentawai
Islands, Indonesia, vol. 441. Geological Society London Special Publications,pp. 47–70. https://doi.org/10.1144/SP441.3, 1.
Hanks, T.C., Kanamori, H., 1979. A moment magnitude scale. J. Geophys. Res. 84 (5),2348–2350. https://doi.org/10.1029/JB084iB05p02348.
Honjo, Y., 2011. Challenges in geotechnical reliability based design. Proceedings of the
3rd International Symp. Geotech. Safety Risk 11–27.
Imamura, F., Yalciner, A.C., Ozyurt, G., 2006. Tsunami modelling manual. http://www.
tsunami.civil.tohoku.ac.jp/hokusai3/E/projects/manual-ver-3.1.pdf.
Ishikawa, Y., Kameda, H., 1988. Hazard-consistent magnitude and distance for extended
seismic risk analysis. In: Proceedings of Ninth World Conference on Earthquake
Engineering. Tokyo-Kyoto, Japan, pp. 89–94.
Japan Hydrographic Association, 2015. M7000 digital bathymetric chart, Japan nauticalchart web shop. http://www.jha.or.jp/shop/.

Great East Japan earthquake. The method was employed to estimate theexceedance probabilities with reference to the tsunami trace heights
observed in a field survey and attained overall consistent results.
Furthermore, we quantified the contribution ratios of different uncertainties and found that the fault slip and modeling error have muchgreater influences than the rake on the overall uncertainty of thetsunami height.
The key feature of the proposed method is the use of an RS obtained
from a series of numerical simulations that are computationally expensive. By performing an MCS with the RS, we can efficiently estimate theprobability density function of the target output. In addition, the proposed method enables us to utilize advanced numerical simulations thatare computationally expensive because these simulations were needed
only to construct the RSs. Furthermore, we can determine engineering
priorities by estimating the contribution of each uncertainty to the target
output. The proposed PTHA method is flexible and allows us to assess
various tsunami risks, for example, not only tsunami heights but also
inundation areas. The novel method will also be useful for assessing thetsunami risks of future earthquakes along other coastlines that haveexperienced tsunami-induced damage in the past.
Future work includes the extension of the proposed method toconsider other uncertainties associated with tsunami hazards. Additionally, it is worth investigating adaptive RS models (e.g., Wong et al.,2005) to use an increasing number of numerical simulation runs toefficiently build a surrogate model of the numerical simulation for
PTHA. Another research direction is to use the proposed approach tovalidate numerical models involving notable geological indicators oftsunami runup and inundation.
Declaration of competing interest
The authors declare that they have no known competing financial
interests or personal relationships that could have appeared to influence
the work reported in this paper.

11

T. Kotani et al.

Coastal Engineering 160 (2020) 103719

Japan Society of Civil Engineers, 2002. Tsunami assessment method for nuclear power
plants in Japan. http://committees.jsce.or.jp/ceofnp/node/.
Japan Society of Civil Engineers, 2011. The method of probabilistic tsunami hazard
analysis (in Japanese).
Kajiura, K., 1963. In: The Leading Wave of a Tsunami, vol. 41. Bulletin of the Earthquake
Research Institute, University of Tokyo, pp. 535–571, 3.
Kajiura, K., 1982. In: Tsunami Energy in Relation to Parameters of the Earthquake Fault
Model, vol. 56. Bulletin of the Earthquake Research Institute, University of Tokyo,
pp. 415–440, 3.
Kanamori, H., 1977. The energy release in great earthquakes. J. Geophys. Res. 82,
2981–2987. https://doi.org/10.1029/JB082i020p02981.
Larsen, B.E., Fuhrman, D.R., 2019. Full-scale CFD simulation of tsunamis. Part 1: model
validation and run-up. Coast. Eng. 151, 22–41. https://doi.org/10.1016/j.
coastaleng.2019.04.012.
Lorito, S., Selva, J., Basili, R., Romano, F., Tiberti, M.M., Piatanesi, A., 2015.
Probabilistic hazard for seismically induced tsunamis: accuracy and feasibility ofinundation maps. Geophys. J. Int. 574–588. https://doi.org/10.1093/gji/ggu408,
200.
McGuire, R.K., 1977. Seismic design spectra and mapping procedures using hazard
analysis based directly on oscillator response. Earthq. Eng. Struct. Dynam. 5 (3),211–234. https://doi.org/10.1002/eqe.4290050302.
McGuire, R.K., 2004. Seismic Hazard and Risk Analysis. Earthquake Engineering
Research Institute, Berkeley.
Mitsoudis, D.A., Flouri, E.T., Chrysoulakis, N., Kamarianakis, Y., Okal, E.A., Synolakis, C.
E., 2012. Tsunami hazard in the southeast aegean sea. Coast. Eng. 60, 136–148.
https://doi.org/10.1016/j.coastaleng.2011.09.004.
Mori, N., Goda, K., Cox, D., 2017. Recent process in Probabilistic Tsunami Hazard
Analysis (PTHA) fo Mega Thrust Subduction Earthquakes. In: The 2011 Japan
Earthquake and Tsunami: Reconstruction and Restoration. Springer, pp. 469–485.
Mori, N., Takahashi, T., 2012. Nationwide post event survey and analysis of the 2011
Tohoku earthquake tsunami. Coast Eng. J. 54, 1250001. https://doi.org/10.1142/
S0578563412500015.
National Research Council, 1988. Probabilistic Seismic Hazard Analysis. The National
Academies Press, Washington, DC. https://doi.org/10.17226/19108.
Necmioglu, O., Ozel, N.M., 2014. An earthquake source sensitivity analysis for tsunami
propagation in the eastern mediterranean, apecial issue on undersea natural hazards.
Oceanography 27 (2), 76–85. https://doi.org/10.5670/oceanog.2014.42.
Okada, Y., 1992. Internal deformation due to shear and tensile faults in a half-space. Bull.
Seismol. Soc. Am. 82 (2), 1018–1040.
Omira, R., Baptista, M.A., Matias, L., 2015. Probabilistic tsunami hazard in the northeast
atlantic from near- and far-field tectonic. Pure Appl. Geophys. 172, 901–920.
https://doi.org/10.1007/s00024-014-0949-.
Omira, R., Matias, L., Baptista, M.A., 2016. Developing an event-tree probabilistic
tsunami inundation model for NE atlantic coasts: application to a case study. Pure
Appl. Geophys. 173, 3775–3794. https://doi.org/10.1007/s00024-016-1367-z.
Park, H., Cox, D.T., 2016. Probabilistic assessment of near-field tsunami hazards:
inundation depth, velocity, momentum flux, arrival time, and duration applied to
Seaside, Oregon. Coast. Eng. 117, 79–96. https://doi.org/10.1016/j.
coastaleng.2016.07.011.
Park, H., Cox, D.T., Barbosa, A.R., 2018. Probabilistic Tsunami Hazard Assessment(PTHA) for Resilience Assessment of a Coastal Community. Nat. Hazards. https://doi.org/10.1007/s11069-018-3460-3.
Park, H., Alam, M.S., Cox, D.T., Barbosa, A.R., Lindt, J.W. van de, 2019. Probabilistic
seismic and tsunami damage analysis (PSTDA) for the Cascadia Subduction Zone

applied to Seaside, Oregon. Int. J. Disaster Risk Reduct. 35, 101076. https://doi.org/
10.1016/j.ijdrr.2019.101076.
Power, W., Downes, G., Stirling, M., 2007. Estimation of tsunami hazard in New Zealand
due to south American earthquakes. Pure Appl. Geophys. 164, 547–564. https://doi.
org/10.1007/s00024-006-0166-3.
Qin, X., Motley, M.R., Marafi, N.A., 2018. Three-dimensional modeling of tsunami forceson coastal communities. Coast. Eng. 140, 43–59. https://doi.org/10.1016/j.
coastaleng.2018.06.008.
Sarfaraz, M., Park, A., 2016. SPH numerical simulation of tsunami wave forces impinged
on bridge superstructures. Coast. Eng. 121, 145–157. https://doi.org/10.1016/j.
coastaleng.2016.12.005.
Sarri, A., Guillas, S., Dias, F., 2012. Statistical emulation of a tsunami model forsensitivity analysis and uncertainty quantification. Nat. Hazards Earth Syst. Sci. 12
https://doi.org/10.5194/nhess-12-2003-2012, 2003–2018, 2012.
Satake, K., Fujii, Y., Harada, T., Namegaya, Y., 2013. Time and Space distribution of
coseismic slip of the 2011 Tohoku earthquake as inferred from tsunami waveform
data. Bull. Seismol. Soc. Am. 103, 1473–1492. https://doi.org/10.1785/
0120120122.
Selva, J., Tonini, R., Tiberti, M.M., Romano, F., Grezio, A., Melini, D., Piatanesi, A.,
Basili, R., Lorito, S., 2016. Quantification of source uncertainties in seismic
probabilistic tsunami hazard analysis (SPTHA), geophys. J. Intell. 205, 1780–1803.
https://doi.org/10.1093/gji/ggw107.
Sraj, I., Mandli, K.T., Knio, O.M., Dawson, C.N., Hoteit, I., 2014. Uncertainty
quantification and inference of Manning’s friction coefficients using DART buoy data
during the To‾hoku tsunami. Ocean Model. 83, 82–97.
Stefanakis, T.S., Contal, E., Vayatis, N., Dias, F., Synolakis, C.E., 2014. Can small islands
protect nearby coasts from tsunamis? An active experimental design approach. In:
Proceedings of the Royal Society A: Mathematical, Physical and Engineering
Sciences, vol. 470. Royal Society publishing, p. 20140575. https://doi.org/10.1098/
rspa.2014.0575.
Suppasri, A., Imamura, F., Koshimura, S., 2012. Probabilistic tsunami hazard analysis
and risk to coastal populations in Thailand. J. Earthq. Tsunami 6 (2), 1250011.
https://doi.org/10.1142/S179343111250011X.
The Port and Airport Research Institute, Nationwide Ocean Wave information network
for Ports and Harbours (NOWPHAS), http://nowphas.mlit.go.jp/nowphasdata/static/sub311.htm.
Thio, H.K., Somerville, P., Ichinose, G., 2007. Probabilistic analysis of strong groundmotion and tsunami hazards in Southeast Asia. J. Earthq. Tsunami 1 (2), 119–137.
https://doi.org/10.1142/S1793431107000080.
Tohoku Earthquake Tsunami Joint Survey Group, 2011. Nationwide field survey of the2011 off the pacific coast of Tohoku earthquake tsunami. J. Jpn. Soc. Civ. Eng. Ser.
B2 (Coast. Eng.) 67, 63–66. https://doi.org/10.2208/kaigan.67.63.
Volpe, M., Lorito, S., Selva, J., Tonini, R., Romano, F., Brizuela, B., 2019. From regional
to local SPTHA: efficient computation of probabilistic tsunami inundation maps
addressing near-field sources. Nat. Hazards Earth Syst. Sci. 19, 455–469. https://doi.
org/10.5194/nhess-19-455-2019.
Wong, S.M., Hobbs, R.E., Onof, C., 2005. An adaptive response surface method forreliability analysis of structures with multiple loading sequences. Struct. Saf. 27 (4),287–308. https://doi.org/10.1016/j.strusafe.2005.02.001.
Working Group on California Earthquake Probabilities, 1995. Seismic hazards insouthern California - probable earthquakes, 1994 to 2024. Bull. Seismol. Soc. Am.
85, 379–439.

12