Journal of Japan Association for Earthquake Engineering, Vol.8, No.4, 2008

DESIGN EARTHQUAKE GROUND MOTIONS FROMPROBABILISTIC RESPONSE SPECTRA: CASESTUDY OF NEPAL
Hari PARAJULI1, Junji KIYONO2 , Yusuke ONO3 and Takahiro TSUTSUMIUCHI4
1

Ph. D. student, Dept. of Urban Management, Kyoto University, Kyoto, Japan,
parajuli@quake2.kuciv.kyoto-u.ac.jp
2
Member of JAEE, Associate Professor, Dept. of Urban Management, Kyoto University, Kyoto, Japan,
kiyono@quake.kuciv.kyoto-u.ac.jp
3
Member of JAEE, Assistant Professor, Dept. of Urban Management, Kyoto University, Kyoto, Japan,
ysk@quake2.kuciv.kyoto-u.ac.jp
4
Master student, Dept. of Urban Management, Kyoto University, Kyoto, Japan,
tsutsumiuchi@quake2.kuciv.kyoto-u.ac.jp

ABSTRACT: Probabilistic hazard estimate through out Nepal considering historicalearthquakes, intra plate slip and faults are done. As a typical case probabilistic spectra areplotted for Pokhara city. For the city, design earthquakes for three probabilities of
exceedences are simulated which can be useful to design new structures and retrofit ofexisting structures.
Key Words: Faults, historical earthquakes, kernel, response spectra, design earthquakes

INTRODUCTION
Nepal takes approximately half length of greater Himalaya, which is part of the trans-alpine belt,
regarded as one the main earthquake prone zones of the world. Because of tectonic movement, many
earthquakes have occurred in this region in its history. Records noted on some Nepalese religious
tracts indicate that a big earthquake hit Kathmandu in June 1255 AD. The quake killed approximately
one-third of its population at that time. Since then severe earthquakes have been reported which
occurred in 1405, 1408, 1681, 1810, 1833, 1866 and 1934 AD (BECA 1993, Ambraseys and Douglas2004). Evidences of a big earthquake in central part of Nepal occurred in the period from 700 to 1100
AD have been published recently but exact occurred date has not conformed yet (Lave et al. 2005).
Historical records show that at least as far back as the early 18th century, damaging earthquakes haveoccurred in the Himalayan region in every few decades. But since 1950, the damaging earthquakes
have not been reported in the region and in some areas. From recent studies depending upon thevarious analyses and GIS (Geographical Information System) data, slip rate of Indian and Tibetan
plates is 19 mm per year (Jouanne et al. 2004). But the calculated slip for whole Himalyan-arc fromoccurred earthquakes is only one third of the observed slip (Bilham and Ambrasseys 2005). It showsthat either the earthquake records are missing or severe earthquakes may be overdue. Many
earthquakes struck greater Himalayan region in past but the worst may be yet to come and it may
occur in or near Nepal. According to the new analysis by Bilham and Ambrssseys 2005, one or more- 16 -

massive earthquakes measuring greater than M8 on Richter scale may be overdue in the Himalaya,threatening the millions of people that live in the region.
More than ten thousand people were killed in 1934 earthquake. Since then, Nepal’s population has
doubled and urban population has increased by a factor of more than ten since last great earthquake. In
urban, most of the people live in poorly constructed houses without considering the seismic codes. In
rural area almost all people live in low strength masonry houses constructed by mainly by stones andbricks which have been proved the major cause of live loss in earthquakes. Considering past human
tolls from earthquakes, population increases that have occurred since then and the added low quality
houses, the future scenario of deaths and damages could be worse than 2005 Pakistan earthquake.
Previous seismic hazard estimate (BECA 1993) is based on uniform distribution of assumed
earthquakes data over big area without well explained details. That study does not consider theearthquakes greater than magnitude 8 whereas their catalogue also has earthquakes greater than thatmagnitude. Thus, revising all earthquake data, selecting suitable attenuation law, peak ground andspectral accelerations are estimated and contours are plotted for the region. For a site, contribution toparticular value of acceleration by various magnitudes and distances are disaggregated. For all bins
significant durations are calculated. Weighted average duration and envelope function are obtained bymultiplying disaggregated hazard by significant duration. Design earthquake is estimated for
calculated weighted average duration fitting with obtained probabilistic spectra. The whole procedureis shown in flow chart (Fig. 1).

Probabilistic seismic hazard analysis

Calculation of probabilistic spectra

Deaggregation of
hazardCalculation of weighted average duration

Simulation of design earthquakes

Definition of envelope
functionFig. 1 Flow chart for simulation process

MAGNITUDE FREQUENCY RELATIONSHIP
Usual practice to develop earthquake frequency relationship is to collect earthquakes, assign to nearest
fault and make the relationship. The site consists of 92 faults (BECA 1993). Only few earthquakes areavailable, most of the faults find very few data and even some of them are empty (Fig. 2). Available
earthquakes are arranged in various magnitudes and years (Table 1). All the faults and historicalearthquakes greater than magnitude M3 around Nepal are plotted in Fig.2. They are scattered andconcentrated out side the faults. Being faults are historical evidences; they can not be neglected eventhough earthquakes have not occurred in time span of data.
Instead of considering individual faults, area sources have been modeled and earthquake density atparticular area cell is calculated considering both historical earthquakes and faults which lay in the cellby kernel method which is explained later. Earthquake catalogue was formed merging the data from
U.S. Geological Survey, National earthquake Information Centre (NEIC), BECA 1993, Ambraseys andDouglas 2004 and Lave et al. 2005. Though historical earthquakes occurred in the region have been

- 17 -

documented since 1100AD, earthquakes data have gone missing (Bilham and Ambrasseys 2005 andJouanne et al. 2004) and available data are not sufficient to justify the current slip rate 19mm/year(Jouanne et al. 2004). Magnitude frequency relationship (eq. 2) has been developed (Parajuli et al.
2007) by calculating slope of the relationship from existing earthquakes and considering 50%
intra-plate slip which clearly satisfies (Fig. 3) the rate of smaller and great earthquakes which is
explained in the following paragraphs.
32
31
30
Latitude

China
29

Nepal
*
Pokhara

28
27
26
25
76

M3-M4
M4.1-M5
M5.1-M6
M6.1-M7
M7.1-M8
M>8
78

India

80

82

84
Longitude

86

88

90

Fig. 2 Historical earthquakes (points) and faults (line)
Table 1 Arrangement of earthquake events
Year
Magnitude
Start End 4.0-4.9 5.0-5.9 6.0-6.9 7.0-7.9 8.0-8.9 Total
1100 1931
1
2
5
5
2
15
1932 1936
1
2
1
1
1
6
1937 1941
2
2
1942 1946
1
1
1947 1951
1
1
1952 1956
4
2
61957 1961
22
1962 1966
7
5
1
13
1967 1971
11
3
1
15
1972 1976
13
3
16
1977 1981
12
12
1982 1986
14
14
1987 1991
29
2
1
32
1992 1996
19
1
1
21
1997 2001
23
2
25
2002 2006
16
16
Total
146
30
12
6
3
197

- 18 -

92

Since the earthquake data have been reported in different magnitudes and intensity scales, all data
were converted to moment magnitude (Hank and Kanamori 1979) using various relationships(McGuire 2004) and scaling relationship for Himalayan region (Ambraseys and Douglas 2004). The
total 197earthquake events with magnitude greater than M4 are grouped as shown in table 1. The data
are not uniformly distributed temporally. Analysis for temporal completeness was performed (Stepp
1972), with events grouped into small intervals of time. If k1, k2, k3, ..kn, are the number of quakes perunit time interval, then an unbiased estimate of the mean rate of earthquakes per unit time interval of
the sample exceeding each magnitude is given by eq. 1. The rates of earthquakes exceeding eachmagnitude using eq. 1 are calculated (Table 2) and plotted in semi log scale (Fig. 3). The magnitude
frequency relationship is represented by eq. 2.

1 n
∑ ki
n i =1
Log ( N / Y ) = 3.47 − 0.76M
rateM =

(1)
(2)

Table 2 Mean rate of historical earthquakes
Magnitude
exceeding
M4
M5
M6
M7
M8

Year
Start
1962
1932
1657
1255
1100

End
2006
2006
2006
2006
2006

Duration(years)
32
75
350
752
602

Average
rate
3.97
0.37
0.069
0.013
0.003

1

Nos.of earthquakes(M>m)/year

10

Log(N/Y)=3.47-0.76M
Log(N/Y)=3.60-0.76M
0

10

-1

10

-2

10

-3

10

3

4

5

6
7
Magnitude
Fig. 3 Magnitude frequency relationship89The slope of magnitude frequency relationship (eq. 2) is less than 1 and shows high rate ofoccurrences of larger magnitude earthquakes. Though historical earthquakes in Nepal have been
documented since 1100 AD, various researchers (Bilham et al. 2007, Feldl and Bilham 2006, Bilham
and Ambrasseys 2005, Jouanne et al. 2004) have pointed out that an intra-plate slip deficit exists in the
Himalayan region. The slip velocity between India and southern Tibet as measured by GPS
measurement is 16-18 mm/year (Bilham et al. 1997) and the slip velocity of central and eastern Nepal
is 19mm/year (Jouanne et al. 2004). The slip rate estimate (Bilham and Ambraseys 2005) based on 500
years of earthquake data accounts for only one third of the total seismic rate of the region, and theremaining slip deficit is equivalent to five earthquakes greater than M8.5.

- 19 -

The question therefore arises of whether this fit is sufficient to represent the seismicity of the region
when a large slip deficit exists and the faults are sufficiently mature to sustain renewed rupture (Feldl
and Bilham 2006) effectively, there is the possibility of a mega thrust event with 600km by 80km
rupture area and more than 9m slip (Bilham and Ambraseys 2005). Slip rate ( s& ) can be obtained frommoment rate ( M& 0 ), shear modulus ( µ ) and fault area ( A f ) using eqs. 3-6 (McGuire 2004).

s& =

M& 0
µAf

(3)

ν kβ exp[β (M min + 16.05 / 1.5)] 1− β / γ
M& 0 = M min
M 0 max − M 01−minβ / γ
γ −β

(

[

k = 1 − e − β (M max − M min )

νM

min

= 10

(a − bM min )

)

]

−1

(4)
(5)
(6)

Where γ = 3.454 and β = 2.303b , and M 0 max and M 0 min are maximum and minimum moments
corresponding to maximum and minimum magnitudes, calculated using eq. 7 (Hank and Kanamori
1979).
log10 M 0 = 1.5M + 16.05(7)
Considering the lower threshold magnitude M5 and maximum magnitude M8.8 (Lave et al. 2005),
fault area 600km by 80km, a and b (coefficients) values from eq. 2 and µ=3.3E11dyne/cm2, the slip
rate is 37% of total slip (1.9cm/year). Alternative explanations are possible for the apparent
contradiction between observed data and slip velocity. Earthquake data may be missing, the magnitude
value of the events may have been underestimated, the slip may be aseismic creep, or the current study
period may not be long enough compared to the recurrence rate of great earthquakes. Efforts havebeen made to find missing data – for example, recently an earthquake in 1100 AD has been identifiedin eastern part of Nepal (Lave et al. 2005) and another one in Garhwal Himalaya (Kumar et al 2006).
Observing eq. 2, slip depends upon a, b values, area of fault and maximum magnitude. The main
shortcoming of this magnitude frequency fit is that it is based on relatively few data. However,
earthquake data are almost complete at magnitude M4, and the period of study since 1100 AD may not
be sufficient for events greater than M8. As three earthquakes greater than M8 have been identified
and the return period for these events is nearly 500 years (Feldl and Bilham 2006), the seismic hazard
assessment should be based on the rates of these events. In that sense, considering the slope (b value)of the magnitude frequency relation and maximum magnitude are constant, and taking the same area
as Bilham and Ambraseys 2005, an increased slip rate to 50% yields a (coefficient of magnitude
frequency relation) value of 3.60, which closely matches the rate of occurrence of events greater than
M8 and complete M4 (Fig. 3, blue line), so that the magnitude frequency relation can be represented
by eq. 8.
Log ( N / Y ) = 3.60 − 0.76M(8)PROBABILISTIC SEISMIC HAZARD ASSESSMENT (PSHA)
Seismic hazard curve for each site can be obtained from summing up the mean rate of exceedences ofall small source zones. Ns is the numbers of sources in the region, total mean rate exceedences(Kramer 1996) for the region is estimated by eq. 9.
Ns Nr Nmν y = ∑∑∑ν i*i =1 j =1 k =1

M min

ρi P[Y > y* | m, r ]P[M = m]P[R = r ]∆m∆r- 20 -(9)ρi = Nλ (m, x )i

(10)

s

∑ λ (m, x )

i

[

i =1

]

where P Y > y * | m, r is conditional probability that ground motion parameter Y (acceleration)exceeds a particular value y* using attenuation law for a given magnitude M and distance R, andP[M = m] and P[R = r ] are probabilities of occurrence of particular magnitude and distancerespectively. No specific attenuation law has been developed for Himalayan region. Atkinson and
Boore 2003 compiled all earthquakes database added many recent earthquakes data from Japanthrough 2001, formed four times bigger database than others’ for subduction zone events anddeveloped new ground motion relation. It is the latest and includes all types of earthquakes and isselected in this study. Earthquake data has magnitude, distance and depths and can be directly used inattenuation law. For area sources, magnitude dependent parameter delta (Atkinson and Boore 2003) istaken as depth. Lengths of faults are taken from BECA 1993. Since, some of the faults do not hold any
earthquake data, maximum magnitude calculated (Wells and Coppersmith 1994) for each faults wereassigned as one equivalent earthquake while calculating density. Here, M and m are used as randomvariable and specific value for magnitude. Total area of 600kmx600km around each site is taken anddivided into 120x120 cells. Distances between centre of cells and site are calculated. Only the cells
within 300km radius are considered in the study. Magnitude is divided into 0.5M and distance into
5km intervals. Nm and Nr are the total numbers of magnitudes and distances bins. Earthquake densities(ρi) for each cell are considered using kernel estimation methods (Woo 1996). The mean activity
rate λ (m, x ) satisfying r ≤ h m j , at a cell x is taken as a kernel estimation sum considering the( )contribution of N events inversely weighted by its effective return period which can be obtained from
eqs. 11-13.
N

λ (m, x )i = ∑j =1

K (m j , rj )T (rj )

 D   h(m j )
K (m, r ) j =  2πh(m j )  rj h(m j ) = H exp(Cm j )

(11)
2− D

(12)
(13)

where, K (m, x ) is kernel function, T (r ) is return period of the event located at distance r, h(m )is kernel band width scaling parameter shorter for smaller magnitude and vice-versa, which is regarded
as fault length, D is fractal dimension, is taken as 1.7 , H and C are constants equals to 1.45 and 0.64.

PROBABILISTIC SPECTRA
Mean rate of exceedences of peak ground accelerations and spectral accelerations at various periods ofvibration from all small sources are obtained from eq. 1 using Atkinson and Boore 2003 attenuation
equation. Assuming earthquake occurrences obey Poisson’s process, peak ground accelerations andspectral accelerations corresponding to 40%, 10% and 5% probabilities of exceedences in 50 years are
determined. Distribution of peak ground accelerations over the country for the three different
probability of exceedences which are equivalent to 98, 475 and 975 years return period respectively
have been plotted in Figs.4-6. Figs. show that eastern and far western regions may experience bigger
peak ground accelerations than in central regions. Two big (greater than M8) earthquake have occurredand lower magnitude historical earthquakes (Fig. 2) are also densely populated in these regions thanin central region. Long faults are also in western and eastern side, however small faults are in central
part. As a typical example, for Pokhara city, latitude 28.2N and longitude 83.98E, peak groundaccelerations and spectral accelerations at various periods for soft soil condition are obtained (Fig. 7).

- 21 -

80
60

1 00
80

80

80

80

220

80

80

100

40

27

60

26.5
80

81

82

100

100 0
8

28 4
0
27.5

60
80

80

*
Pokhara60

Latitude

80

60

60

28.5 60

0
12

100

0
12

0
10
100

120

30.51
4
2108
30 000 60
129.5 121040160
10
0
29

83

84
85
Longitude

60

80

86

87

88

Fig. 4 Peak ground acceleration (gal) in 40% in 50 years (soft soil)

200

45
0

250

27
26.5
80

35
0

500

500

*
Pokhara

81

82

0
45

400
350

350

83

400

450
400

27.5

350

50
25 330
0 0

5 50

300

40
0

Latitude

28

50
40 5
400
0

500

400

400

0
45

28.5

0
45

0
40
350 30029

0 0
50 4 5

0 0
40 35

600

29.5 500550

550

600

550
66500
30 0

0
55

30.5

84
85
Longitude

86

87

88

Fig.5 Peak ground acceleration (gal) in 10% in 50 years (soft soil)

00
10
929.5 800
00
29

80
0

600

0
60

0
70

500

* Pokhara
700

70
0

60
0

700

800

800
70
0

400

27

700

27.5

800

600

60
0
28 500

700

28.5

70
0

0
70

Latitude

30

90
0

1000

0 0
90 80

0
80

900

0
90

30.5

600

0
60

26.5
80

81

82

83

84
85
Longitude

86

87

88

Fig. 6 Peak ground acceleration (gal) in 5% in 50 years (soft soil)- 22 -

1200

40% in 50 yrs

Acceleration (gal)1000

10% in 50 yrs
5% in 50 yrs

800
600
400
200
0
0.0

0.5

1.0

1.5

2.0

2.5

3.0

Natural Period (secs.)

Fig. 7 Response spectra (5% damping) at Pokhara (soft soil)

DURATION AND ENVELOPE FUNCTION
Duration of earthquake is function of magnitude and epicentral distance, thus each earthquake has
separate duration. A probabilistic spectrum consists of many earthquakes. Even within duration span
acceleration amplitudes are not uniform. Thus total duration (TD) is divided into three parts as shownin Fig.8. In the first part up to TB acceleration will be ascending, between TB to TC, it is much effective
and after TC, it starts descending. TB and TC are calculated using equations following Osaki 1994 and
J

1

K

E(t)
0.1
OTB

B

TC C

TD

D

0.1

Contibution (%)

Contribution (%)

Fig. 8 Division of duration and envelope function

0.05

0
0

1 M(5.0-5.5)
2 M(5.6-6.0)
3 M(6.1-6.5)
4 M(6.6-7.0)
5 M(7.1-7.5)
6 M(7.6-8.0)
7 M(8.1-8.5)
8 M(8.6-9.0)

0.2

0.1

0
0
100
Distance (km)200

100
Distance (km) 200
300

Magnitude

40% in 50 years

300 Magnitude
10% in 50 years

Fig. 9 Disaggregation PGA at Pokhara (soft soil)- 23 -

TD is calculated using equation following Kemption and Stewart 2006. It is significant duration which
is defined as the time interval across which 5 to 95% of total energy is dissipated.

TB = [0.12 − 0.04(M − 7 )]TD
Tc = [0.50 − 0.04(M − 7 )]TD

(14)
(15)

  exp(b + b (M − 6 ))  −13

1
2



101.5 M +16.05
ln TD = ln
+ c2 r + c1s 
4.9.106 β





(16)

where, b1, b2, c1, c2 are coefficients equal to 2.79, 0.82, 1.91, 0.15 respectively, β is shear wave
velocity equal to 3.2km/sec. s is soil type and equal to 1 for soil, and zero for rock, M is magnitude ofearthquake and r is epicentral distance. Now, envelope function E(t) is calculated using eqs. 17-19.

0 ≤ t ≤ TB : E (t ) =  t 
 TB 
TB ≤ t ≤ Tc : E (t ) = 1

2

(17)
(18)

 ln 0.1
(t − TC )
0 ≤ t ≤ TB : E (t ) = exp TD − TC(19)

In order to find out appropriate duration, all accelerations are disaggregated in terms of magnitudes
and distances (Fig. 9). Disaggregated cells have separate magnitudes, distances and weights. Separate
significant durations for all cells using eqs. 14-16 were calculated, then, obtained durations are
multiplied by corresponding cell’s weight obtained from disaggregation and weighted durations arefound summing up all bins’ durations. We obtain duration for particular value of response acceleration
corresponds to specific period. However, we need a single value applicable for all acceleration that fallin specified probability of exceedence in specified period of years. Thus, using similar procedure,
durations for all accelerations are calculated. These calculated durations have combined effects of
distance, magnitude and weight of deaggregation. For lower accelerations, contribution of lower
magnitude earthquakes in hazard is significant but duration is short. For distant earthquakes duration islong but its contribution to hazard is low. For frequent earthquakes, smaller magnitude earthquake also
contribute more but when probability decreases contribution of higher magnitude earthquake increases.
In an average, the duration remains almost same for all earthquakes which fall in same return period.
Thus, the weighted durations corresponding to same return period are in close margin. Then weighted
average duration is calculated from all accelerations. For spectra shown in Fig. 7, durations arepresented in Table 3.
Table 3 Weighted average duration for Pokhara
DurationsTD
TB
TC

40% in 50 yrs
30.0
3.5
15.0

10% in 50 yrs
45.0
5.3
23.0

5% in 50 yrs
51.0
5.9
26.0

DESIGN EARTHQUAKES
Design earthquakes are simulated for calculated durations (Table 3). Total duration is first divided into
small interval dividing by number N as shown in eq. 20. Using incremental time, Fourier transformpair Ck and xm are evaluated through eqs. 21-24.

- 24 -

TD
N

(20)

1 N −1∑ xm exp − i 2πkm N , k=0, 1, 2, … , N-1
N m =0

(21)

∆t =

( (
))
x = ∑ C exp(i (2πkm )), m=0, 1, 2, … , N-1
N

Ck=

N −1

m

k =0

k

22)

Ck = Fk (cos φk + i sin φk )

(23)

Rd =

S ds
S dt

(24)

where Fk is Fourier amplitude,φ is phase angle, Rd, is ratio of simulated (Sds) to target spectra (Sdt).
Envelope plot (Fig.8) was divided into small increments. Using eqs. (17-19) value of envelope
function E(t) was calculated. Cumulative value of E(t) was calculated at each interval and these allvalues are divided by final sum which give cumulative probability density function from zero to one.
Considering total duration as three hundred sixty degree, phase angles for every time steps were
determined randomly using probability density function. Phase angles can be obtained from previous
earthquake records. However, there are no recorded earthquake acceleration histories for the region.
Thus, random phase angles estimated from envelope function were used. At first, Fourier amplitude Fk
is assumed unity. From, Fast Fourier Transform (FFT), accelerations at each interval was calculated.
Accelerations were obtained from Inverse Fourier Transform. To make the simulated ground motion
similar to natural earthquakes, the acceleration obtained from inverse Fourier transform were
multiplied by envelope function. Using the calculated accelerations, response spectra was determined
and compared with original spectra called target spectra. Ratio from simulated spectra to target spectra
at each interval was obtained. New Fourier amplitude was then calculated multiplying old amplitude
by obtained ratio. Again, FFT were calculated and accelerations were determined. The process was
repeated until the simulated spectra and target spectra fit well. The schematic diagram of wholeprocess is shown in Fig. 10.

Durations TB, TC, TD

Envelope function and phase angles
Fk=1,Rd=1
Fourier Amplitudes (Fk=RdFk)

Fourier transform (Fk, x(t))

Plotted spectrafits target spectra

No

Rd =

Yes

S ds
S dt

Design earthquakes
Fig. 10 Flow chart for generating acceleration histories- 25 -

100
50
0
-50
-100
0

10 15 20 25 30
Duration (secs.)a) 40% in 50 yrs. (random wave input)

5

750
500
250
0
-250
-500
-750
0

Acceleration (gal)Acceleration (gal)

450
300
150
0
-150
-300
-450
0 5 10 15 20 25 30 35 40 45
Duration (secs.)b) 10% in 50 yrs.(random wave input)

Acceleration (gal)Acceleration (gal)The accelerations histories obtained from this procedure for three probabilities of exceedences havebeen plotted in Fig. 11. From simulated acceleration histories acceleration response spectra were
calculated and plotted against target spectra (Fig.12). For target peak ground accelerations are taken atzero period. But, for simulated spectra, plotting starts from 0.1 sec. for 975 and 475 years return
periods and 0.05 sec. for 98 years return period. However, the spectra from simulated earthquake andoriginally calculated from attenuation law are in good agreement. In real earthquakes, accelerations
start from zero, increase gradually, attains peak values and decrease to zero finally. Simulated
earthquakes also look similar to natural earthquakes. Amplitudes and nature of simulated groundmotions are totally depend on target spectra.

10

20
30
40
50
Duration (secs.)c) 5% in 50 yrs. (random wave input)

750
500
250
0
-250
-500
-750
0

10 20 30 40 50
Duration (secs.)d) 5% in 50 yrs. (Kobe 1995 input)

Fig. 11 Simulated acceleration history for various return periods

Spectral acceleration (gal)1500

1 Target (RT 98 yrs.)
2 Target (RT 475 yrs)
3 Target (RT 975 yrs.)
4 Simulated (RT 98 yrs)
5 Simulated (RT 475 yrs.)
6 Simulated (RT 975)
7 Simulated (RT 975 yrs.)

1000

500

0
0

0.5

1

1.5
2
2.5
Period (secs.)
Fig. 12 Comparison of simulated and target spectra

3

In order to see the differences of simulated ground motion from real earthquake and randomly
occurred earthquake, a simulated earthquake ground motion for 5% probability of exceedence in 50- 26 -

years has been shown in Fig. 11-d. Input was given Kobe 1995 earthquake. When it is fitted with 5%
probability of exceedences in 50 years, the original earthquake completely changes its shapes
according target spectrum. Thus, whatever is the input either real earthquake random occurrences there
are many possibilities, however, the maximum amplitude of acceleration in both case are almost
similar (Fig.11 c-d).

DISCUSSION
Probabilistic seismic hazard estimate for a region where earthquakes data are missing has been done.
There are many small faults around the big Himalayan thrust. But sufficient historical data are notavailable and available earthquake data does not satisfy the current intra-plate slip. Many of the faults
are empty of historical earthquakes. However, they are the surface evidences of rupture occurred in its
history. The method employed here combines three parts; historical earthquakes, faults and intra plate
slip and estimate the probable hazard for the region. The density of earthquakes (Fig. 2) and contour
plot (Fig.4-6) shows that the eastern and western part may experience larger accelerations than thecentral part. Then probabilistic spectra for the city Pokhara has been plotted (Fig. 7). Probabilistic
earthquakes were obtained for weighted average duration which obtained from significant duration
multiplied by weights of deaggregation and taken average for all durations corresponding to specified
return period. Previous researchers have also given efforts to estimate risk consistent ground motion
for example Kameda and Nojima 1988. Simulation of ground motion has been carried out focusing thecharacteristic or significant earthquakes. The method employed here through deaggregation is different
than others. If we look at the deaggregation plot (Fig. 9), we can get magnitude and distance forsignificant earthquake corresponding to specified return period as a whole. The method employed here
considers through the entire possible earthquakes that fall in the specified highest weight in the plot.
Even though it is the highest, it may have just 10% or less weightage in overall hazard. Thus the bins
other than highest weight have significant contributions. Few significant earthquakes lack to represent
probabilistic earthquake for return period and simulate equivalent earthquake. Thus it differs fromother estimates which are based on significant or characteristic earthquakes. This method can be usefulto estimate hazard where sufficient data are not available.

CONCLUSION
Peak ground acceleration for three probabilities of exceedences were calculated and contours havebeen plotted. Accelerations are far higher than previously calculated by BECA 1993. Thus Revision of
hazard estimate is necessary to incorporate in seismic codes. As a typical example probabilistic spectra
for Pokhara city have been plotted. For the city, three separate probabilistic earthquakes were
simulated. These earthquakes can be used for dynamic analyses for different life span structures. In
usual practice, acceleration histories obtained from previous earthquakes are taken. However, they can
not represent certain probability of exceedences. Same earthquake can not be used for varying life
span structures. Thus, either new structures or existing structures, they are designed for specific life
span. Thus, these probabilistic earthquakes can be useful to design for new structures or retrofitting ofexisting structures.
REFERENCES
Ambraseys, N.N. and Douglas, J. (2004). “Magnitude calibration of north Indian earthquakes”,
Geophysics. J. Int., Vol. 159, 165-206.
Atkinson, G.M. and Boore, D.M. (2003). “Empirical ground-motion relations for subduction-zone
earthquakes and their application to Cascadia and other regions”, Bull. Seismol. Soc. Am., Vol. 93(4), 1703-1729.

- 27 -

Bilham, R., and Ambraseys, N. (2005). “Apparent Himalayan slip deficit from the summation ofseismic moments for Himalayan earthquakes”, 1500-2000, Current science, Vol. 88 (10),1658-1663.
BECA World International (New Zealand) in association with SILT Consultants (P.) Ltd. (Nepal),
TAEC Consult (P.) Ltd. (Nepal), Golder Associates (Canada) and Urban Regional Research (USA)(1993). “Seismic Hazard Mapping and Risk Assessment for Nepal”.
Feldl, N., and Bilham, R. (2006). “Great Himalayan earthquakes and the Tibetan plateau.” Nature
Publishing Group, 444(9), 165-170.
Hank, T.C., and Kanamori, H. (1979). “A moment magnitude scale.” Journal of Geophysics Res., 84,2348-2350.
Jouanne, F., Mugnier J. L., Gamond, J.F., Le Fort, P., Pandey, M.R., Bollinger, L., Flouzat, M.,
Avounac, J.P. (2004). “Current Shortening across the Himalayas of Nepal”, Geophysics. J. Int., Vol.
157, 1-14.
Kameda, H. and Nojima, N., (1988). “Simulation of risk-consistent earthquake motion”, Earthquake
Engineering and Structural Dynamics, Vol. 16, 1007-1019.
Kempton, J. J. and Stewart, P. J. (2006). “Prediction equations for significant duration of earthquake
ground motions considering site and near-source effects”, Earthquake spectra, Vol. 22(4),985-1013.
Kramer, S.L. (1996). “Geotechnical earthquake engineering”, Prentice-Hall International series in
Civil Engineering and Engineering Mechanics.
Kumar, S., Wesnousky, S.G., Rockwell, K.T., Briggs, R.W., Thakur, V.C, and Jayangondapernumal R.
(2006). “Paleoseismic evidence of great surface-rupture earthquakes along the Indian Himalaya.”
Journal of Geophysical research, 111.
Lave, J., Yule, D., Sapkota, S., Basnet, K., Madan, C., Attal, M., and Pandey, R. (2005). “Evidence fora great medieval earthquake (1100 A.D.) in the central Himalayas”, Science www.sciencemag.org,
Vol. 307.
Mcguire, R. K. (2004). “Seismic hazard and risk analysis.” Earthquake Engineering Research Institute,
MNO-10Parajuli, H., Kiyono, J., Scawthorn, C. (2007). Probabilistic seismic hazard assessment of Kathmandu,
Proceedings of twentieth KKCNN symposium on civil engineering, Jeju, Korea
Stepp, J.C. (1972). “Analysis of completeness of the earthquake sample in the Pudet Sound area andits effect on statistical estimates of earth
Wells, D.L. and Copersmith, K.J. (1994). “New empirical relationships among magnitude, rupture
length, rupture width, rupture area, and surface displacement.” Bull. Seismol. Soc. Am., 84, (4),974-1002.
Woo, W. (1996). “Kernel estimation methods for seismic hazard area source modeling”, Bull. Seismol.
Soc. Am., Vol. 86 (2), 353-362.
大崎順彦 1994. 新・地震動のスペクトル解析入門 鹿島出版会(Submitted: May 30, 2008)(Accepted: September 18, 2008)- 28 -