NUMO-TR-08-03

Development of Methodologies
for the Identification of Volcanic
and Tectonic Hazards to Potential
HLW Repository Sites in Japan- The Tohoku Case Study -Neil Chapman1, Mick Apted2, John Beavan3, Kelvin Berryman3,
Mark Cloos4, Charles Connor5, Laura Connor5, Olivier Jaquet6,
Nicola Litchfield3, Sue Mahony7, Warwick Smith3, Steve Sparks7,Mark Stirling3, Laura Wallace3
1

MCM Consulting, Switzerland
Monitor Scientific LLC, USA
3
GNS Sciences, New Zealand4
University of Texas, USA
5
University of South Florida, USA
6
In2Earth Modelling Ltd, Switzerland
7
University of Bristol, UK
2

March 2009
Nuclear Waste Management Organization of Japan (NUMO)

2009 年 3 月

初版発行

本資料の全部または一部を複写・複製・転載する場合は，下記へ
お問い合わせください。
原子力発電環境整備機構 技術部〒108-0014 東京都港区芝 4 丁目 1 番地 23 号 三田 NN ビル 2 階
電話 03-6371-4004（技術部） FAX 03-6371-4102
Inquiries about copyright and reproduction should be addressed to:
Science and Technology Department
Nuclear Waste Management Organization of Japan (NUMO)
Mita NN Building 2F, 4-1-23, Shiba, Minato-ku, Tokyo 108-0014, Japan
Phone +81-3-6371-4075 FAX +81-3-6371-4083©原子力発電環境整備機構(Nuclear Waste Management Organization of Japan) 2009

NUMO-TR-08-03

Development of Methodologies
for the Identification of Volcanic
and Tectonic Hazards to Potential
HLW Repository Sites in Japan- The Tohoku Case Study -Neil Chapman1, Mick Apted2, John Beavan3, Kelvin Berryman3,
Mark Cloos4, Charles Connor5, Laura Connor5, Olivier Jaquet6,
Nicola Litchfield3, Sue Mahony7, Warwick Smith3, Steve Sparks7,Mark Stirling3, Laura Wallace3
1

MCM Consulting, Switzerland
Monitor Scientific LLC, USA
3
GNS Sciences, New Zealand4
University of Texas, USA
5
University of South Florida, USA
6
In2Earth Modelling Ltd, Switzerland
7
University of Bristol, UK
2

March 2009
Nuclear Waste Management Organization of Japan (NUMO)

Contents
1

Introduction .......................................................................................................................... 1

2

The Geology and Tectonics of the Tohoku Region............................................................. 4

3

2.1

Northern Honshu – Tohoku .......................................................................................... 6

2.2

Basement of Tohoku .................................................................................................... 6

2.2.1

Forearc region ....................................................................................................... 6

2.2.2

Volcanic arc region ................................................................................................ 7

2.2.3

Backarc region ...................................................................................................... 7

2.3

Active Deformation: Earthquakes and Faulting ........................................................... 7

2.4

Example Locations ....................................................................................................... 8

Evaluation of Rock Deformation in Tohoku ....................................................................... 10
3.1

Strain Rates from Surface Deformation ..................................................................... 11

3.1.1

Active fault source map and logic tree ................................................................ 11

3.1.2

Uplift and tilt data and development of logic trees .............................................. 13

3.1.3

Active fault strain map and example location calculations .................................. 15

3.1.4

Tilt strain map and example location calculations ............................................... 18

3.1.5 Combined Surface Deformation Strain Hazard ...................................................... 21
3.2

3.2.1

What exactly is GPS measuring?........................................................................ 24

3.2.2

Extracting the component of the GPS signal relevant to rock deformation hazards
26

3.2.3

Application of our proposed method to the Tohoku case study .......................... 27

3.2.4

GPS velocities, uncertainties and logic tree construction ................................... 29

3.2.5

One representative elastic block model for northern Honshu ............................. 31

3.2.6

GPS residual upper plate strain map and example location calculations ........... 33

3.3

4

Strain Rates from Global Positioning System (GPS) Data ........................................ 24

Strain Rates from Seismicity Data ............................................................................. 38

3.3.1

Seismicity Dataset ............................................................................................... 39

3.3.2

Logic Tree Structure ............................................................................................ 39

3.3.3

Results and Discussion ....................................................................................... 41

3.4

Strain Rate Comparisons ........................................................................................... 45

3.5

From Strain Rate to Hazard ....................................................................................... 47

3.5.1

Hazard Calculation .............................................................................................. 47

3.5.2

Probability of Exceedance of Threshold Strain Rates......................................... 51

Volcanic Hazard Analysis .................................................................................................. 57
4.1

Introduction to the methodology ................................................................................. 57

4.1.1

Geological Models ............................................................................................... 58

4.1.2

Data Uncertainty.................................................................................................. 60

4.1.3

Probabilistic Assessment of Volcanic Hazard ..................................................... 61

4.2

Alternative Data Sets and Cladistics .......................................................................... 62i

4.2.1

Method ................................................................................................................ 63

4.2.2

Results ................................................................................................................ 63

4.2.3

Clade-volcano cluster correlations ...................................................................... 63

4.3

4.3.1

Spatial Density .................................................................................................... 65

4.3.2

Estimating Spatial Density with Kernel Methods ................................................. 68

4.3.3

Cox process with multivariate potential ............................................................... 68

4.4

5

6

Probability Models ...................................................................................................... 65

Application of Probability Models in Tohoku .............................................................. 70

4.4.1

Kernel Density Estimation ................................................................................... 71

4.4.2

Cox Process Model ............................................................................................. 72

4.5

Assessment of Probability of Volcanic Disruption at 14 Example Locations ............. 74

4.6

Sensitivity to Probability Models ................................................................................. 76

Conclusions and Recommendations................................................................................. 81
5.1

Rock Deformation....................................................................................................... 81

5.2

Volcanism ................................................................................................................... 81

5.3

Closing Remarks ........................................................................................................ 82

References ........................................................................................................................ 84

Appendix
The Fourteen Example Locationsii

1 Introduction
NUMO is responsible for the siting, development and operation of a deep geologicalrepository for high level waste (HLW) in Japan. The process is expected to take at least 15
years to reach the point of repository construction. During the period before this, NUMO willneed to evaluate sites that emerge from the ‘volunteer process’ (whereby local communities
have been invited to volunteer to be considered as potential hosts for the repository) andselect a preferred site. This evaluation will involve initially surface based and thenunderground site characterisation work. Underground characterisation work will only takeplace at the preferred site. Prior to the surface based investigations, volunteer sites will havefirst had to pass a test of general suitability and NUMO will then have carried out a detailed,
literature-based preliminary evaluation of suitability, prior to accepting them as ‘Preliminary
Investigation Areas’ (PIAs). Because Japan lies in such a tectonically active region of theworld on the Pacific rim (the so called ‘ring of fire’), a key aspect of all these steps isconsideration of the susceptibility of a site to future tectonic activity and tectonically driven
processes and events.
In particular, the potential for volcanic and rock deformation impacts on a repository siteneeds to be considered at each stage of NUMO’s siting programme. Whilst the nationwide
evaluation factors for qualification (EFQs) for PIA acceptance are designed to remove clearly
unsuitable sites from consideration, they cannot guarantee that, over the next tens of
thousands of years, the risks of tectonic hazard for a chosen PIA will be acceptable. This isbecause large parts of Japan that are potentially suitable for siting are directly affected tovarying extents by rock deformation, the peripheral impacts of volcanic activity or the
possibility of new magma intrusion or volcanic activity. The EFQs were only intended as a‘blunt instrument’ to prevent obviously poor candidates entering the siting process.
Consequently, additional ‘sharper’ and more refined techniques are required to evaluate sitesthat pass the EFQ test, so that NUMO can have a clear idea of the likelihood and potential
impacts of tectonic events and processes at each PIA. The ITM project1 has been designed toprovide NUMO with such a methodology, based upon state-of-the-art approaches usedinternationally, developed and extended for the specific purposes of NUMO and the specificconditions of Japan: hereafter, we refer to it as the ‘ITM Methodology’.
The ITM methodology is essentially probabilistic in nature. A probabilistic approach is seen by
the ITM expert group as the only realistic means of addressing the uncertainties in predicting
possible hazards when there is marked variability in the spatial distribution, the timing, theintensity and the style of the volcanic and deformational events and processes being
evaluated. The probabilistic approach being developed is based upon and strongly supported
by deterministic models of the underlying tectonic processes that lead to magma intrusion,
volcanism and rock deformation.
NUMO is developing both the probabilistic ITM methodology and other, independent,
deterministic approaches, in parallel projects which will eventually be deployed at volunteer
sites when they arise. The weight that will be given to deterministic and probabilistic
evaluation results will depend, to some extent, on the nature and the geographical location ofthese sites.
The probabilistic ITM methodology is able to be used at three important stages of NUMO’ssiting programme:•

STAGE 1: during the literature survey (LS) stage when potential PIAs are being
assessed. The ITM methodology will use currently available information to allow

1

The ITM project was started by taking a recommendation in the second International Tectonics Meeting in 2004 that
NUMO should prepare a probabilistic methodology as well as deterministic ones which will mainly be employed insite evaluation.

1

comparison of sites in terms of confidence that they are likely to prove acceptablewith respect to tectonic impacts.
•

STAGE 2: during the planning of the PIA site investigations, to identify geoscientificinformation requirements that will be needed to refine the Stage 1 analysis.

•

STAGE 3: at the point where PIAs are being evaluated and compared in order toselect a preferred site (or sites) for detailed investigation (as DIAs).

The ITM project is mainly concerned with Stages 1 and 2. Application of the methodology in
Siting Stage 3 is several years into the future and it is expected that it will be most efficient tocarry out any necessary updates/refinements on a region-specific basis during the PIA
investigations when NUMO has narrowed down to a group of sites. The ITM project involves
methodology development and testing only and does not include actual deployment forvolunteer sites/regions.
The overall structure of the ITM methodology is described in Chapman et al. (2009 in press)and consists of:•assembling nationally available data and alternative models of the nature, causes andlocations of tectonic processes and events;•using probabilistic techniques to evaluate the likelihood and scale of future tectonic
processes and events, shown as a function of their type and geographicaldistribution;•

feeding information on these potential likelihoods and impacts to NUMO’s
performance assessment team so that feedback can be provided on repository
performance under tectonic stress;•

providing clearly justified and traceable input to decision-making on consequent site
suitability.

For convenience, the methodology for rock deformation and volcanic hazards assessment
has been applied as two parallel tasks. This recognises the fact that, although the concept ofeach approach as shown above is similar, in some parts of the methodology they differ
significantly in detail. Consequently, it was found that two teams with different specialities(structural, geophysics and tectonics specialists; volcanologists) worked efficiently in parallel.
However, it is most important that, if carried out in this manner in future, the two ‘discipline’
teams integrate their work frequently (more so than has been possible in this development
project) as there are clear overlaps in the processes being evaluated (e.g. magma intrusion
has an impact on rock stress regimes and vice versa). NUMO will need to ensure that suchintegration is carried out effectively when the methodology is applied to ‘real’ sites.
The broad structure presented above is shown in more detail in the top-level methodology‘road-map’ in Figure 1.1. It comprises a series of eight Steps, distinguishing in the early Steps
between the ‘rock deformation’ evaluation and the ‘volcanic’ evaluation, where they involve
significantly different activities.
This Report covers the final results of a Case Study designed both to help develop and thento test the ITM Methodology. The Case Study was carried out using information about the
Tohoku region of northern Honshu. This area was selected for the study as it is the centre ofcurrent discussions about the varied strain response of the crustal plate to subduction of the
Pacific Ocean plate (the key current tectonic driver for much of Japan) and about themechanisms that underlie the apparent clustering of Quaternary volcanoes in much of
Honshu. A second Case Study is adapting the ITM methodology developed in the Tohoku
area for examining the very different tectonic situation of Western Japan – Kyushu, and theresults of this work will be reported separately.

2

Select region of interest
that includes one ormore potential
repository sites(~104 – 105 km2)

Data & Conceptual Models
Existing catalogues
and databases

Step 1
Assemble regional Volcanic,
Deformation & Tectonic Data

Step 2• evaluate data, filter, obtain ground truth• identify alternative conceptual models
for magma intrusion & rock deformation

Rock Deformation: Step 3b
Strain Models

Magma Intrusion: Step 3aClassify & Group

Cladistic
Analysis

Logic
Trees

 strain models for each of several
independent sets of data into classes of event

ExpertElicitation

 produce alternative groupings

 in each dataset, weight eachalternative model whencalculating strain
integrated
analysis

 time series analysis of each

Step 4a

Step 4b

Intrusion Probability Maps

Probabilistic Strain Maps
• for each independent dataset• difference maps

for each Group and time set

Step 5Integrated Evaluation foreach site

Step 6: Interface with PA• identify impact scenarios for PA• assign probabilities to each

Step 7
SitingDecision
MakingAssign ‘Siting Confidence’ Levels
•
•
•
•
•

Repository
Concepts

PA evaluations(outside ITM project)for each magma intrusion Groupfor each type of rock deformation
low tectonic susceptibility area
uncertain suitability for some RCs
potentially unsuitable for any RC

Step 8
PI fieldplanning

Define PIA and DIA field datarequirements• verify Step 6• reduce uncertainties

Figure 1.1: The Steps in the ITM methodology, shown as a top-level road-map.

3

2 The Geology and Tectonics of the Tohoku Region
Japan is part of the "Ring of Fire," the belt of earthquakes and volcanic activity that
distinguishes the active margins of the Pacific Ocean from the passive margins of the Atlantic
Ocean. In the 1920s, the great seismologist K. Wadati, discovered that earthquakes beneath
northern Japan form an inclined zone extending from locations very near the Japan Trench todepths of about 500 km beneath the Japan Sea. This trench is a north-south trending
bathymetric depression about 150 km east of the mainland that is as deep as 9000 m. The
towering volcanoes of northern Japan are centred about 75 km from one another and form acurving line that is about 100 km above the inclined seismic zone. In the 1930s, Japaneseseismologists discovered that many of the earthquakes near the trench are the result of large
thrust fault movements indicating the floor of the Pacific basin is moving beneath northern
Japan.
The cause of the active tectonic movements affecting the Earth was not well explained until
the late 1960s, when the realization was made that the outer part of the Earth is divided into
pieces, known as plates, about 100 km thick (Takeuchi et al., 1970). The plates consist ofboth the crust and cool uppermost mantle that has sufficient long-term strength thatearthquake-generating ruptures can occur. Plate boundaries come in three forms that arerecognized by the common type of fault movement: divergent with normal faulting,convergent with thrust faulting, and transform with strike-slip faulting. Plates separate fromone another in the process of seafloor spreading to form ocean ridges and come together inthe process of subduction to form the ocean trenches, inclined seismic zones and curvinglines of explosive volcanoes. They slide past one another at transform boundaries.
Most of the geology of Japan is a result of subduction-related processes since the Mesozoic(Sugimura and Uyeda, 1973). The patterns of active faulting and seismicity along with direct
GPS measurements delineate where tectonic motions are underway (Sagiya et al., 2000).
To a first order, the current tectonics of the Japanese islands can be explained by theinteraction of four plates: Pacific, Philippine, Eurasian and North American (see Figure 2.1).
The eastern part of the Eurasian plate is broken with a large fragment, the Amur subplate,currently moving as a distinct kinematic entity (Wei and Seno, 1998; Heki et al., 1999). The
North American plate continues across the Bering Sea into eastern Asia and down past theKamchatka-Kurile trench segments to Japan. An elongate southern prong of the North
American plate extends southwards to Japan. This prong has broken off and is currently
moving as a distinct kinematic entity, the Okhotsk subplate.
In detail, the current tectonics for southern Kyushu and the Ryukyu arc-trench system is amanifestation of the interaction of the Philippine and Eurasian plates. The current tectonics ofSouthwest Japan (Kyushu, Shikoku, and southwest Honshu) are a manifestation of theinteractions of the Philippine plate and Amur subplate. The current tectonics of Northeast
Japan (northern Honshu and Hokkaido) are a manifestation of the interactions between the
Amur and Okhotsk subplates with the Pacific plate. Subduction along the Japan Trench at a
speed of about 9 cm/year (90 km/Ma) is concurrent with convergence near the eastern edgeof the Sea of Japan at a speed of about 1 to 1.5 cm/year (10 to 15 km/Ma) (Tamaki and
Honza, 1985; Okamura et al.,1995).
Four arc segments merge to form the Japanese islands. The northern half of Honshu island
is a subduction segment commonly referred to as northern or northeast Japan. Hokkaido is
the southern end of the Kurlie trench and arc. East-central Honshu (Izu peninsula) is thenorthern end of the Izu-Bonin-Mariana trench and arc. These three arc-trench segments are
all the product of the westward subduction of the Pacific plate. Kyushu island is the northernend of the Ryukyu trench and arc. The active tectonics of Kyushu and Shikoku islands are
movements related to the northwestward subduction of the Philippine plate. East-central
Honshu, just south of Tokyo, near 34°N is presently the convergent junction of three major
plates.

4

Figure 2.1: The principal plate tectonic features of Japan, showing the location of the Tohoku region ofnorthern Honshu that forms the focus of this Case Study. The location of active faults and volcanoes
that have been active in the last 10,000 years is also shown.

Prior to Eocene time, the entire length of present day Honshu was underthrust by the Pacificplate. At about 45 Ma, tectonic movements changed along the western margin of the Pacific
Ocean. Whether the Philippine plate is a fragment broken off the edge of the Asian or the
Pacific plates is debated. It is unique amongst the large plates in that it is nearly entirely
surrounded by subduction zones. The eastern edge of the Philippine plate is underthrust by
the Pacific plate forming the Izu-Bonin-Mariana arc and trench. The western edge of thePhilippine plate subducts beneath Asia creating the Ryukyu arc and trench. The Izu arc andtrench intersection with the Japan trench has progressively migrated along southwest Japan
to its present position during the late Cenozoic. In so doing, subduction beneath Shikoku
island and southwest Honshu at the Nankai Trough has changed from the underthrusting of
the Pacific plate to the underthrusting of the Philippine plate.

5

2.1

Northern Honshu – Tohoku

Since the advent of plate tectonic theory, Northern Japan, from 35°N to 43°N, has been
widely recognized as the type example of an active subduction zone. The northeast Japansubduction zone, like all others, can be subdivided into the trench axis, trench slope, forearc
basin, volcanic arc, and backarc regions. The tectonic activity associated with subduction
comes in the form of faulting, magmatic intrusions and extrusions and vertical movementsthat range from folding to tilting of fault blocks to regional warpings.
The earthquakes are the result of episodic (stick-slip) movements along faults. Arc volcanismis the result of water-fluxed melting deep in the mantle because subduction injects water. The
melting generates hydrous basalts with elevated concentrations of most of the trace elementsthat are abundant in sediments. These magmas rise buoyantly towards the surfaces andcentres of arc volcanism are located where the magmas work their way through the crust.

2.2

Basement of Tohoku

The basement of Japan is largely composed of Palaeozoic and Mesozoic subduction-relatedrock terranes created during westwards subduction along the eastern edge of Asia. The
northernmost part of Honshu is underlain by the Oshima Belt, a Jurassic accretionary prism.
Central Honshu is underlain by the Ashio Belt, largely composed of Triassic to Jurassic rocks
that are variably deformed and little metamorphosed accretionary prism materials and the
Abukuma Belt which is thoroughly metamorphosed under high-temperature / low-pressure
conditions. Between these two Mesozoic terranes is the Kitakami Belt, a varied mixture of
sedimentary, igneous and metamorphic rocks of Palaeozoic age. This juxtaposition makes it
appear that the early Mesozoic history of subduction and accretionary prism growth was
followed by an episode of transform faulting (Taira, 2001). Strike-slip faulting shuffled theaccretionary prism - forearc basin terranes (Abukuma-Ashio and Oshima belts) with the arcbasement terrane (Kitakami Belt) that had formed along the edge of the Asian continent.
The present phase of westward dipping subduction began at about 130 Ma (Engebretson etal., 1985). Cretaceous and Cenozoic intrusives and volcanics of the northern Honshu arc are
emplaced into and on top of these basement terranes. Early Cretaceous plutonics are mostabundant near the eastern coastline and occur offshore (Finn, 1994). Late Cretaceous arc
magmatism appears to have been concentrated near the west coast of Honshu. Thesemagmatic rocks were down-dropped during the rifting associated with the opening of the Sea
of Japan and deeply buried by younger sediments (Finn et al., 1994). The volcanic arc
appears to have established it present position since about 20 Ma. Crustal thickness beneathmost of northeast Honshu is between 30 to 35 km.
Most of the plutonic rocks were emplaced at depths less than 10 km as the feeder roots of
volcanoes very similar and perhaps identical to those dominating the landscape of today.
Most ancient volcanoes were destroyed by erosion long ago, but relicts of various ages arefound in many locations. In most of northern Japan, the plutonic rocks are emplaced into
older plutons and volcanics or a basement of sedimentary rocks that was deformed andmetamorphosed to varying degrees.

2.2.1

Forearc region

The forearc region is directly underlain by the subducting Pacific plate. The geology of theforearc region of northeast Japan records a long history of subsidence and volcanogenic
sediment accumulation. Because the upwards flow of heat through the descending plate
beneath the forearc block is slower than the speed of subduction, geothermal gradients arevery low across the forearc region (Honda, 1985). High-pressure/low-temperature blueschist
facies metamorphism is occurring at depth. The crystalline basement underlying the forearc
region is cold and strong. Consequently the sediments in the overlying forearc basin are little
deformed.

6

The forearc region of northeast Japan has a rather smooth slope down to the Japan Trench.
The Cenozoic tectonic history of the forearc region has largely been one of non-accretion andslow subsidence resulting from subduction erosion of the base of the hanging wall block
above the descending Pacific plate (von Huene and Lallemand, 1990; Heki, 2004).

2.2.2

Volcanic arc region

The area around active volcanic arcs is a region of high heat flow because heat is advected toshallow levels in magmas. Recent seismic tomography studies indicate most of the active
volcanic centres in northeast Japan are centred above areas with slow seismic velocities in
the upper mantle and crust (Tamaru et al., 2002). The current volcanic arc position wasestablished at about 20 Ma (Kondo et al., 2004). Since the mid-Miocene, numerous
magmatic centres have formed and decayed. The lifespan of typical large arc-type composite
volcanoes appears to be about 1 to 2 million years (Davidson and De Silva, 2000). Volcanicarcs are noted for giant composite volcanoes. But clusters of satellite cones and flows are
typical. Along the northeast Japan arc, ten volcanic clusters are about 50 km wide that are
separated by gaps between 30 to 75 km wide (Tamura et al., 2002).
Many subduction zones, including northern Japan, have scattered occurrences of volcanismin the backarc. Whether these occurrences are a manifestation of limited magma generationin the underlying mantle or limited pathways for intrusion to the surface is a matter of debate.
In some cases, as is currently occurring in the Okinawa Trough behind the Ryukyu Arc, thebackarc is a region of extensive volcanism because of seafloor spreading. Where this occurs,
the amount of volcanism in the backarc region will greatly exceed the activity along the arc.

2.2.3

Backarc region

The geology of the backarc and arc region of northeast Japan records a history of profoundchange. Until mid-Cenozoic time, the crust forming the bulk of the basement of Japan was
the edge of Asia. The Japan Sea formed by continental rifting which evolved into seafloor
spreading with ocean crust formation between 23 to 14 Ma (Jolivet et al., 1994). This rifting
event was centred along the volcanic arc that had formed a line delineated by the westernedge of Honshu.
Volcanic arcs can become centres of divergence that evolve into sea floor spreading centres.
This has occurred at several sites along the western margin of the Pacific basin. Major
episodes of backarc spreading occurred along Izu-Bonin-Mariana and Tonga-Kermadec
subduction zones. In the backarc regions of these oceanic subduction zones, lines of extinct
volcanoes parallel the current volcanic arc with a region of new ocean crust in between.
It is evident that rifting is localized to the part of the lithosphere that was sufficiently weakened
by the ascent of magmas, so that a change in the force balance near the subduction zone led
to lithospheric divergence that is centred on the arc. Backarc spreading is a subordinate
process to subduction. During the time seafloor spreading created a 200 km width of new
ocean crust beneath the Japan Sea, on the order of 2000 km width of Pacific plate wasconsumed by subduction at the Japan Trench.

2.3

Active Deformation: Earthquakes and Faulting

Since 1900, nineteen earthquakes of M 7 to M 8.2 have occurred directly beneath the coast to
about 50 km offshore northeast Japan between depths of about 15 to 60 km (Kawakatsu andSeno, 1983). These events occur along the plate interface as a direct manifestation of theunderthrusting of the Pacific plate. The updip cut-off is known as the seismic front with little
seismicity occurring seaward in the hanging wall block. The origin of the abrupt seismic front
along Northeast Japan and elsewhere is debated. The cessation of downdip seismicity is
generally considered a manifestation of where the material along the plate interface issufficiently heated so that it flows rather than fractures.

7Less frequently, earthquakes occur beneath the land. Most occur at depths less than 15 km(Zhao et al., 2000). Seismicity is sparse and shallow in areas of high heat flow and nearly
absent beneath centres of active volcanism.
A distinct belt of seismicity and surface faulting extends along the western coast of northeast
Japan north and south of Niigata. Scattered reverse-slip earthquakes as large as M7.8
nucleate at depths of 10 to 20 km (Okubo and Matsunaga, 1994). Faulting and folding in thebelt is a result of slow convergence between Amur and Okhotsk subplates. Most of thisconvergence is localized along the western edge of northern Honshu. This area of active
folding, reverse faulting, uplift and erosion was an area of subsidence during the Miocene
phase of rifting.
Normal faults created during rifting to form the Japan Sea were buried by highly volcanogenic
sediment shed westward into the widening depression. The convergent deformation in the
Niigata region began in the Pliocene, at about 3 Ma (Okamura et al., 1995). The crystalline
basement in this region primarily deforms by reversing the movement on normal faults
created during the rifting phase – ‘inversion tectonics’. The overlying sediments respond to
the deformation primarily by drape folding over the rising fault blocks. Geometricalcomplexities and new breakouts occur in many areas because the old fault system does noteverywhere have optimal orientations for the imposed movements.
Some convergence is occurring within the region of the volcanic arc. These movements havecreated the Backbone Range, a topographic divide that has profoundly affected regional
drainage and sedimentation patterns (Sato, 1994). Between volcanic centres, some of theshortening occurs by episodic earthquake-generating stick-slip along reverse faults (Zhao etal., 2000). Directly beneath most volcanic centres, the temperatures are so high that the crustis ductile and there is little seismicity. The uplift of the Backbone Range began since about 2
Ma (Sato and Amano, 1991).
In short, the tectonic environment of northern Japan is dominated by subduction at the Japan
Trench causing arc magmatism with second-order convergent motions that are concentrated
in the backarc (Niigata) region. The current mode of tectonism in northern Honshu began atabout 4 Ma and has been well established since 2 Ma (Taira, 2001). The initiation of thistectonic regime is probably a manifestation of the creation and movement of the Amur andOkhotsk subplates. Why these convergent movements began, which are secondary
compared to fast, long-term subduction at the Japan Trench, is not known but they are part ofa modest change of plate motions around the Pacific basin that occurred between 5 to 3 Ma
(Cox and Engebretson, 1986; Pollitz, 1986).

2.4

Example Locations

For the purposes of the present study, a group of fourteen locations (‘example locations’) was
selected across the Tohoku region so that the results of the probabilistic rock deformation andvolcanic evaluations could be discussed in terms of different tectonic settings within theregion. Within the ITM methodology, the two evaluations produce quantitative estimates of
future upper crustal strain and the likelihood of new magmatic intrusion, discretised into 5 x 5km blocks across the whole region. Identifying example locations allows these quantitative
estimates to be presented and discussed for varying tectonic regimes in the Tohoku region.
The example locations selected are shown in Figure 2.2. The Appendix to this report provides
a brief description to each location, as well as showing the location-specific quantitative
results of the deformation and volcanism evaluations. The evaluations themselves are
discussed in detail in the following Sections (Sections 3 and 4).

8Figure 2.2: The fourteen example locations used in the Tohoku Case Study (green dots). Also shown
are volcanic edifices (blue triangles) and active faults (bold, black lines). See Appendix for details.

93 Evaluation of Rock Deformation in Tohoku
The ITM approach to developing a methodology for rock deformation hazard analysis in the
Tohoku case study has required the integration and comparison of several datasets, each ofwhich contains information on past rock deformation and, therefore, potential future hazard.
These datasets are:(i) a 400 year long seismicity record;(ii) a 10-100 thousand year long record of past large earthquake contained in active fault
studies, and deformation captured in a surface deformation map compiled from upliftand folding of fluvial and marine terraces;(iii) a c.10 year record of high precision, horizontal, surface deformation captured bycontinuous GPS velocity data, that is consistent with global plate motion models
averaged over the past several million years.
To investigate multiple lines of evidence for rock deformation a common measure between
the datasets is required. The rate of strain in the crust of the earth provides that common
measure. Strain rate describes deformation (e.g. 10-6/year is the annual change in length of a
1 km long section of the earth’s surface by 1 mm). The strain rate can be used to estimate,through probabilistic methods and assumptions regarding the relationship between strain rateand earthquake occurrence, the expected strength of earthquake shaking, amount of tectonic
displacement, or other tectonic hazard that could directly impact a repository over time
periods in the order of 102 – 106 years.
Each of the strain rate datasets has relative strengths and weaknesses. The surface
deformation data provide estimates of strain for geological time periods, but is thought to be
an incomplete dataset with respect to Class C faults (<0.1 mm/year slip rate). The uplift andderived folding deformation comes from a single, preliminary compilation by Tajikara (2005).
The geodetic model has excellent spatial coverage but represents a very short time period of
observation, and is very sensitive to modelling parameters (e.g. the degree of interseismic
coupling assumed for the subduction zone of the Japan Trench). The historical seismicity
model utilises the c. 400 year historical record of Japan, but there are considerable
uncertainties associated with the recurrence parameters calculated for this model. The
following section provides a brief description of the datasets and methods used to derivestrain rates from three datasets. Well-established seismological and geodetic methods havebeen used to convert the three suites of data into strain rate.
The methodology that is developed is probabilistic. Probabilistic hazard assessment has been
the standard for the quantification of seismic, and to a lesser extent other tectonic hazards forthe nuclear industry for over 20 years, and there is now a standard set of guidelines fordevelopment of a PSHA for a nuclear facility (SSHAC, 1997). We have introduced
probabilistic methods to the High Level Waste repository siting project for two principal
reasons. Firstly, it provides a consistent basis for comparison between sites or regions(particularly pertinent at the literature evaluation stage for NUMO), and secondly, aprobabilistic approach facilitates a robust assessment of uncertainty. The most likely scenario
events for a site can be determined from the probabilistic models by disaggregation of thehazard curve to find the largest contributors to the hazard. A key component of probabilistic
hazard assessments is the inclusion of numerous viable tectonic/conceptual models. An
expert elicitation process is commonly used to determine which of a suite of alternativemodels should be used in the development of the logic tree. The rock deformation groupconducted an expert elicitation workshop in Tokyo in August 2005, and obtained input from apanel of Japanese experts in GPS (Manabu Hashimoto, Kyoto University; Satoshi Miura,
Tohoku University; Takuya Nishimura, Geographical Survey Institute), in seismology (Naoshi
Hirata, ERI, Tokyo; Kunihiko Shimazaki, ERI, Tokyo), and active faulting/surface deformation(Takashi Nakata, Hiroshima Technical University; and Yoko Ota, Emeritus Professor of
Yokohama National University; Yuichi Sugiyama, Director of Active Fault Research Centre,

10Geological survey of Japan). The panel helped us to develop and weight branches of a logic
tree for each of the rock deformation datasets.
Two types of uncertainty are incorporated into probabilistic modelling. Data commonly containsome statistical components of uncertainty such as the “a” & “b” values of the earthquake
frequency-magnitude distribution or fault slip rates where detailed studies have been
completed. These statistical uncertainties are referred to as aleatory uncertainty. Any
parameter can therefore be described by the median of a log-normal distribution. The
distribution may be characterised by standard deviations. To account for knowledge
uncertainty an analysis must incorporate the full spectrum of alternative credible scientificmodels (e.g., different models and/or values of model parameters) that would be provided if
the informed scientific community could be convened to provide their input directly.
Furthermore, since the analysis should determine the distribution of credible estimates, it is
implied that a measure of the scientific credibility of individual assessments shall beconsidered. This is accomplished by weighting the results based on their scientific credibility.
In the following sections we discuss the data available for the case study, the methods wehave applied to the data and the strain models and hazard assessment we have developed.
We present the results in both map form for the whole of the study area, and also for 14
example locations that serve to illustrate the rock deformation and volcano methodologies atpotential site scale.

3.1

Strain Rates from Surface Deformation

Surface fault traces, uplifted and subsided former shorelines, and deformed formerly flat or
inclined surfaces such as river terraces, are all manifestations of past rock deformation, often
associated with major earthquakes resulting from fault rupture in the brittle shallow crust.
Fault slip rates are readily converted to strain rates provided the width over which thedeformation occurs (the fault’s process zone) can be determined. Uplift or subsidence ofmarine terraces often, but not always, results from movement on active faults in the coastal
area. By assuming fault dip then uplift can be converted to fault slip rate and thus to strain
rate. Deformation of marine and fluvial terraces and other areal measures of deformationsuch as elevations of former erosion surfaces are readily assessed as surface strain rates
provided the age of the deformed feature can be ascertained.
In this project we develop an integrated surface deformation strain model by combining active
fault data with folding, uplift and subsidence, as defined by tilted fluvial and marine terraces.

3.1.1

Active fault source map and logic tree

The onshore active fault sources (Figure 3.1) were taken from the active fault databaseprovided on the website by the Active Fault Research Centre (AFRC) of National Institute ofAdvanced Industrial Science and Technology (AIST). The faults in the AFRC database
comprise the 100 most active faults in Japan, so do not include slow slip faults (e.g., Class C
Faults, which have slip rates of 0.01-0.1 mm/year) 2. Some onshore Class C Faults wereadded from data supplied by Prof. Yoko Ota (Figure 3.1), but are considered likely tounderestimate the true number. Offshore, the fault sources were obtained from the summary
map of Japan Sea faults (Ohtake et al., 2002).
The fault source locations were simplified by digitising the centreline of mapped fault traces,
and thus make no assumptions about fault complexity. Instead, as described below, fault
width is now to be treated in the active fault branch of the logic tree (Figure 3.2). In somecases the fault source lengths were extended to take into account possible underestimation of
fault lengths based on fault trace mapping alone (i.e., the preservation of fault traces depends
2

At the beginning of this study in 2004. The present database contains 307 active faults some of which have slip
rates less than 0.1 mm/year.

11

on the preservation of suitable geomorphic surfaces). The extensions are defined along
topographic range-fronts, and generally follow bedrock faults (Nakata and Imaizumi, 2002).
Vertical slip rates for onshore faults were obtained from the AFRC database, where available.
In many cases a range of values are given, but in some cases, only a preferred value,calculated from displacement divided by recurrence interval, is available. For these,
uncertainties still await definition. For other onshore faults, slip rates were calculated frommeasured offsets of Nakata and Imaizumi (2002). For a few onshore faults no data wasavailable, and for these, slip rates were assumed to be the same as other nearby faults forwhich there is some information. Class C faults were assigned the full range of Class C Fault
slip rate of 0.01-0.1 mm/year, and the mean taken as the midpoint (0.05 mm/year). Offshore,
fault slip rates were obtained from Chapter 7 (by Y. Okamura) in Ohtake et al. (2002), andwhere no data are available, slip rates were also assumed to be the same as other nearby
faults, for which there is some information.

Figure 3.1: Active Fault source map for the Tohoku region. Fault sources are colour coded into groups
relating to their data source, as summarised below the map. Vertical slip rates, in mm/year, are labelled
and their data sources are also summarised below the map. The locations of the 14 example locations
are shown in yellow.

Fault width options are defined for two models, a sedimentary rock model (2-20 km), and ahard rock model (0.01-1 km). Two models of dip were also added, which encompass therange of expert opinions of fault dip in the Tohoku area. These are model A (30º-60º) andmodel B (15º-70º).
The order of the branches of the active fault branch of the logic tree is:•for a given fault with predefined segments and lengths•

strain is calculated for a range of widths12•and dips•from vertical slip rate, which is converted to horizontal shortening rate

Weightings are given to each of the branches as assigned by the expert panel (Figure 3.2).

Figure 3.2: Summary logic tree for calculating strain rate from surface deformation. Note only one set ofbranches is shown. The weightings (wt) for the active fault model are as assigned by the expert panel,
but the uplift mechanism model and the tilting model were assigned by the report authors.

3.1.2

Uplift and tilt data and development of logic trees

Surface tilt rates have been derived from the surface uplift contour map for the onshore
Tohoku region since 125 ka (Tajikara, 2005) from pairs of river terraces and from marine

13

terraces (Figure 3.3). Surface uplift has been converted to strain rates from regional tilting
rates, by simply dividing the tilt (mm/km) by the age (125,000 years), noting again that 10-6
per year strain rate = 1 mm/km/year. Given that most of the active faults are approximately
north-south striking and are essentially pure reverse faults, a W-E tilt direction is given thehighest weighting; consistent with the predominant E-W shortening in the Tohoku region.
Thus the main uncertainties, captured by the branches of the logic tree are the tilt (a functionof the uplift measurements) and the age (Figure 3.2). To estimate the surface strain rates
from the 125 kyr uplift data, we fit a surface to the uplift data using the GMT software(http://gmt.soest.hawaii.edu/), using the “surface” function, which is an adjustable tension
continuous curvature surface gridding algorithm. The smoothness of the fit can be adjusted bychanging the “tension factor” used, which a number between 0 and 1 with a tension factor
value of 0 producing the minimum curvature solution. The tilt can then be calculated by taking
the gradient of the surface (here, we used the GMT function, gradient). We considered twocases (see Fig. 3.2 for weightings): (1) the estimated tilt is the maximum gradient at a given
point, and (2) the tilt is the gradient of the uplifted surface in an east-west direction. We also
considered a variety of values for the tension factor (see Fig. 3.2).

Figure 3.3: Uplift (and subsidence) contour map constructed from 125 ka fluvial and marine terraces.
From Tajikara (2005).

14

3.1.3

Active fault strain map and example location calculations

In Figure 3.4 we show an example active fault strain map for the highest weighted realisationof the active fault branch of the logic tree shown in Figure 3.2 (fault dip of 50º, fault processzone width 10 km, and mean vertical slip rate). Active fault strains are calculated for only
those grid squares which contain active faults, which results in very high strain rates in thosesquares, no strain outside of them, and a blocky appearance to the strain map. Alternate
strain maps derived from differing parameters would vary considerably, particularly withestimation of fault width (encompassing fault complexity and the width of the process zone
around an active fault). We also draw attention to the probable incomplete mapping of Class
C faults in the study area, and therefore the likely underestimation of active fault strains.

Figure 3.4: Strain map calculated for the active faults shown in Figure 3.1. The map is for the highest
weighted realisation of the active fault branch of the logic tree shown in Figure 3.2

In Figure 3.5 we show the relationship between known active faults and the 14 example
locations where we will compare rock deformation strain and hazard. The grid is the 5x5 km
we have used in the study, and the effect of fault width choices can be clearly seen onwhether an example location contains active fault strain or not. Seven example locations fallwithin the process zones up to 20 km distant from the delineated fault trace. The strain values
for all of the branches of the active fault branch of the logic tree are shown as a histogram in
Figure 3.6. Location 8 is situated close to one of the most active faults mapped along theJapan Sea coastline (Figure 3.1), and as a result the histogram contains the highest strain
rates. Strain rates on active faults commonly exceed 100 nanostrains/year (100 x 10-9 = 10-7
strain/year), which as we discuss later is indicative of hazard suggesting low confidence orprospectivity of such a location for a HLW underground repository. Locations on or very closeto active faults are already excluded from consideration as a repository under NUMO
regulation, but we include these examples to demonstrate the application of the rockdeformation method in a wide range of settings.

15

Figure 3.5: Locations of the 14 example locations with respect to active faults. The centreline of eachactive fault is shown as a black line, and the shaded buffer zones around each fault indicates the 5km,
10km and 20km wide process zones around each fault.

16

1
0.9

1
0.9

0.8
0.7
0.6

0.8
0.7

0.5
0.4
0.3
0.2
0.1
0
24 47

0.6
0.5
0.4
0.3
0.2
0.1
0

70 93 116 139 162 185 208 231 254 277 300 323

1

Location 3 - FAULTS

Location 13 - FAULTS
1
0.9

0.8
0.7

0.8
0.7

0.6
0.5
0.4
0.3
0.2
0.1
0
24 47

0.6
0.5
0.4
0.3
0.2
0.1
0

70 93 116 139 162 185 208 231 254 277 300 323

170 93 116 139 162 185 208 231 254 277 300 323
nanos trains /ye ar

Location 7 - FAULTS

Location 14 - FAULTS

0.8
0.7

0.8
0.7
0.6

Weighted Frequency

1
0.9

0.6
0.5
0.4
0.3
0.2
0.1
0
24 47

24 47

nanos trains /ye ar

1
0.9

1

70 93 116 139 162 185 208 231 254 277 300 323
nanos trains /ye ar

1
0.9

1

24 47

nanos trains /ye ar

Weighted Frequency

Weighted Frequency

1

Weighted Frequency

Location 11 - FAULTS

Weighted Frequency

Weighted Frequency

Location 2 - FAULTS

0.5
0.4
0.3
0.2
0.1
0
1

70 93 116 139 162 185 208 231 254 277 300 323

24 47

70 93 116 139 162 185 208 231 254 277 300 323
nanos trains /ye arnanos trains /ye ar

Location 8 - FAULTS

Weighted Frequency

1
0.9
0.8
0.7
0.6
0.5
0.4
0.3
0.2
0.1
0
1

24 47

70 93 116 139 162 185 208 231 254 277 300 323
nanos trains /ye ar

Figure 3.6: Histograms of active fault strain at seven example locations located on or close to active
faults, as shown in Figure 3.5. The discrete probability histograms show the probability of a strain rate
being within a 1 nanostrain bin width. The histograms are constructed by sampling all possiblecombinations of parameter values in the active fault branch of the logic tree (Figure 3.2) and calculatinga strain rate for each combination. The histogram is produced by sorting the calculated strain rates
according to 1 nanostrain bins. The horizontal extent of the curve on the x-axis and shape of thehistogram therefore fully quantifies the underlying distribution of strain rate uncertainty. Note that thehistogram for location 8 is truncated, and the total range is 1045 nanostrains/year.

173.1.4

Tilt strain map and example location calculations

As discussed above, we expect that areas on or very close to known active faults are alreadyexcluded as potential sites for a HLW repository. Therefore, data that may provide insights
into the location of slow-slipping faults that have not been mapped, or “hidden” faults, perhaps
buried beneath young alluvium, or locations where there might be “blind” faults, or wherefaults may extend to in the future, are of great interest. Therefore, we have developed thesurface deformation dataset. At present this contains active faults, as discussed in sections
3.1.1 and 3.1.3, and an uplift contour map (Figure 3.3) described in section 3.1.2. The tilt
strain calculations from the uplift contour map dataset are presented in this section and the
combined surface deformation calculations in the next (section 3.1.5).
The tilt strain map derived from the highest weighted path of the logic tree (W-E tilt, tension of
0.15, 125 ka age) is shown in Figure 3.7. High strain rates correlate with close contours of the
uplift map, coinciding with rapid tilt of marine and river terraces. There is also a close
correlation between active fault locations and rapid tilt rates (compare Figures 3.4 and 3.8),
suggesting that known active faults are the driver for some of the surface deformation strain
field. There are other areas where high rates of tilt do not coincide with known active faults,
possibly suggesting hidden, low slip rate, or evolving fault systems. It is also possible that the
uplift dataset is influenced by volcanic-driven deformation which needs considerations,especially at a site level. Areas of low surface tilt frequently correspond to uplifted plateau
areas across Tohoku.

Figure 3.7: Strain map constructed for tilt calculated from the uplift and subsidence contour map in
Figure 3.3. The map shows, for each grid square, the weighted average of all logic tree combinations
shown in Figure 3.2. There are no tilt data offshore - the strain values in the offshore area are the resultof edge effects of curve fitting to the tilt surface - see section 3.1.2.

Two aspects of the strain rates derived from surface uplift data require some cautious
application. Firstly, it is probable that uplift associated with volcanism has been accounted foras tectonic deformation in the model. This aspect of volcano-tectonic interaction in Tohoku
would be a fruitful aspect of further research and probabilistic hazard method development.
Secondly, some uplift with rather broad wavelength may be driven by deep-seated processesof subduction, such as underplating of sediment, and not therefore an indicator of upper plate
fault activity. The tilt dataset is interpreted as a proxy for low slip rate, blind, or hidden faults.
A refinement of the methodology could identify tilt inflection points to indicate where the18localisation of strain is likely to occur (possibly using elastic dislocation theory of Okada, 1971
to suggest the range of possible fault geometry) or evaluate fold wavelength as a criterion forevaluating what the causation process of folding may be. These enhancements of themethodology would be advantageous when additional and increased resolution uplift and tilt
datasets were available for inclusion in the hazard evaluation.
In addition to the map representation of surface tilt strain we have also developed strain rate
histograms for each of the 14 example locations to illustrate the range of strain rate derived
via sampling of the logic tree (Figure 3.8). A characteristic of these histograms is the narrow
range of strain rates calculated. The narrow range of the strain rates results from the very
simple (almost certainly over-simple) formulation of the logic tree. With incorporation ofalternate datasets we would expect the strain rates to attain a broader distribution. At present
it appears there are marked differences between locations indicating a good diagnostic tool
for site comparison, but we need to be cautious because of the limitations noted above. Strain
rates across the various example locations range from almost zero strain to more than 50
nanostrains/year.
Location 1 - TILT

Location 4 - TILT

0.4

0.4

0.35

0.35
0.3
Weighted Frequency

0.25
0.2
0.15
0.1
0.05

0.25
0.2
0.15
0.1
0.05

0
nanostrains/year

89

85

81

77

73

69

65

61

57

53

49

45

41

37

33

29

25

Location 5 - TILT

0.4

0.8

0.35

0.7

0.3

0.6
Weighted Frequency

Weighted Frequency

21

nanostrains/year

Location 2 - TILT

0.25
0.2
0.15
0.1
0.05

0.5
0.4
0.3
0.2
0.1

0
1

5

0

9 13 17 21 25 29 33 37 41 45 49 53 57 61 65 69 73 77 81 85 89

15

9 13 17 21 25 29 33 37 41 45 49 53 57 61 65 69 73 77 81 85 89

nanostrains/year

nanostrains/yr

Location 3 - TILT

Location 6 - TILT

0.4

0.4

0.35

0.35

0.3

0.3
Weighted Frequency

Weighted Frequency

17

9

0

9 13 17 21 25 29 33 37 41 45 49 53 57 61 65 69 73 77 81 85 89

13

5

1

1

5

Weighted Frequency

0.3

0.25
0.2
0.15
0.1
0.05

0.25
0.2
0.15
0.1
0.05

0
89

85

81

77

73

69

65

61

57

53

49

45

41

37

33

29

25

21

17

9

13

5

1

0
1

nanostrains/ye ar59

13 17 21 25 29 33 37 41 45 49 53 57 61 65 69 73 77 81 85 89
nanotrains/year

Figure 3.8 (continued on next page): Tilt strain histograms at 14 example locations from surface
deformation. The discrete probability histograms show the probability of a strain rate being within a 1
nanostrain bin width. See Figure 3.1 for locations, and caption of Figure 3.6 for full explanation of
construction of histograms.

19

Location 8 - TILT

Location 11 - TILT

0.4

0.4

0.35

0.35
0.3
Weighted Frequency

0.25
0.2
0.15
0.1
0.05

0.25
0.2
0.15
0.1
0.05

0

57

61

65

69

73

77

81

85

89

61

65

69

73

77

81

85

89

53

49

45

Location 12 - TILT

0.4

0.4

0.35

0.35

0.3

0.3
Weighted Frequency

0.25
0.2
0.15
0.1
0.05

0.25
0.2
0.15
0.1
0.05

0

53

49

45

41

37

33

29

25

21

17

13 17 21 25 29 33 37 41 45 49 53 57 61 65 69 73 77 81 85 89

13

9

1

5

9

0
1

5

Weighted Frequency

41

nanostrains/yea r

Location 9 - TILT

nanotrains/ye arna notrains/year

Location 10 - TILT

Location 13 - TILT

0.4

0.4

0.35

0.35
0.3
Weighted Frequency

0.3
0.25
0.2
0.15
0.1

0.25
0.2
0.15
0.1
0.05

0.05

0

0
1

5

9

1

13 17 21 25 29 33 37 41 45 49 53 57 61 65 69 73 77 81 85 89

59

13 17 21 25 29 33 37 41 45 49 53 57 61 65 69 73 77 81 85 89
nanostrains/year

nanostrains/yea r

Location 14 - TILT
0.4
0.35
0.3
Weighted Frequency

Weighted Frequency

57

na notrains/ye ar

37

33

29

25

21

17

9

13

1

89

85

81

77

73

69

65

61

57

53

49

45

41

37

33

29

25

21

17

9

13

5

1

0
5

Weighted Frequency

0.3

0.25
0.2
0.15
0.1
0.05
0
1

5

9 13 17 21 25 29 33 37 41 45 49 53 57 61 65 69 73 77 81 85 89
nanostrains/yea r

Figure 3.8: continued from previous page.

20

3.1.5 Combined Surface Deformation Strain Hazard
To calculate a surface deformation strain rate, the combination of active faults and regional
uplift needs to be evaluated. In this study the maximum rate of strain, derived either from
active fault sources (including the width of the process zone) or from the regional tilt model
away from active faults, has been determined. Future development of the surface deformation
dataset could include other marine and river terrace compilations, uplift rate derived from a
range of dating techniques, or long term average rates derived from tilt and elevation of strata
such as shallow water limestone rocks that were originally formed horizontally and near tosea level. The surface deformation dataset is currently restricted to on-land parts of the study
area, but this could be expanded to the offshore using recent basin sedimentary markers andhigh resolution seismic stratigraphy to define relative uplift subsidence contours. Further
stages of expert elicitation would also be desirable to assess the confidence of thosealternate datasets.
Figure 3.9 shows the highest weighted surface deformation strain map. As noted above, thisis the strain rate plotted in each grid square as either the highest weighted branch for active
fault strain, or the highest weighted branch for tilt strain, depending on whichever is greater.
Comparison with Figures 3.4 (the active fault strain map) and 3.7 (the tilt strain map) showsthat the surface deformation strain rates are dominated by the active faults. That is, the map
shows narrow high strain zones along the active faults with lower strain rates in between.
Note that there is no tilt strain data offshore, so zero strain is shown between the active faults.

Figure 3.9: Surface deformation strain map. The strain values have been constructed from two datasets,
active faults (onshore and offshore) shown in Figure 3.1, and tilt calculated from the uplift andsubsidence contour map shown in Figure 3.3. The map shows the weighted average of all realisationsof the logic tree shown in Figure 3.2.

21

In Figure 3.10 the strain values for all branches of the logic tree in Figure 3.2 are shown for
the 14 example locations. Not surprisingly, the plots are generally similar to those for theactive faults (Figure 3.6) or tilt (Figure 3.8). That is, locations which are close to active faults
are generally dominated by active fault strain (which is usually higher than the tilt strain), andvice versa. There are however, some locations which are close to active faults, but where thetilt strain is higher (e.g., 11 and 14). As noted in section 3.1.4, this may indicate the presenceof nearby hidden, low slip rate, or evolving fault systems, but this interpretation should betreated with caution given the initial nature of the tilt dataset. Nevertheless, theseobservations demonstrate the value of combining the datasets to better constrain the surface
deformation strain and as a diagnostic tool for location comparison.

1
0.9

1
0.9
0.8
0.7
0.6
0.5
0.4
0.3
0.2
0.1
0

0.8
0.7
0.6
0.5
0.4
0.3
0.2
0.1
0

24 47 70 93 116 139 162 185 208 231 254 277 300 323

1

Location 2 - SURF DEF

Location 6 - SURF DEF
1
0.9

0.8
0.7
0.6

0.8
0.7

0.5
0.4
0.3
0.2
0.1
0
24 47

0.6
0.5
0.4
0.3
0.2
0.1
0

70 93 116 139 162 185 208 231 254 277 300 323

170 93 116 139 162 185 208 231 254 277 300 323
nanos trains /ye ar

Location 3 - SURF DEF

Location 7 - SURF DEF

0.8
0.7

0.8
0.7

Weighted Frequency

1
0.9

0.6
0.5
0.4
0.3
0.2
0.1
0
24 47

24 47

nanos trains /ye ar

1
0.9

1

70 93 116 139 162 185 208 231 254 277 300 323
nanos trains /ye ar

1
0.9

1

24 47

nanos trains /ye ar

Weighted Frequency

Weighted Frequency

1

Weighted Frequency

Location 5 - SURF DEF

Weighted Frequency

Weighted Frequency

Location 1 - SURF DEF

70 93 116 139 162 185 208 231 254 277 300 323

0.6
0.5
0.4
0.3
0.2
0.1
0
1

nanos trains /ye ar

24 47

70 93 116 139 162 185 208 231 254 277 300 323
nanos trains /ye ar

Figure 3.10 (continued on next page): Strain histograms at 14 example locations from surface
deformation. See Figure 3.1 for locations, and caption of Figure 3.6 for full explanation of construction ofhistograms. Note there are no tilt data at the offshore location 7.

22

1
0.9

1
0.9

0.8
0.7

0.8
0.7

0.6
0.5
0.4
0.3
0.2
0.1
0

0.6
0.5
0.4
0.3
0.2
0.1
0

70 93 116 139 162 185 208 231 254 277 300 323

1

Location 9 - SURF DEF

Location 12 - SURF DEF

0.8
0.7

0.8
0.7

0.6
0.5
0.4
0.3
0.2
0.1
0

0.6
0.5
0.4
0.3
0.2
0.1
0

70 93 116 139 162 185 208 231 254 277 300 323

1

Location 10 - SURF DEF

Location 13 - SURF DEF

0.8
0.7
0.6

Weighted Frequency

0.8
0.7
0.6
0.5
0.4
0.3
0.2
0.1
0

0.5
0.4
0.3
0.2
0.1
0
1

70 93 116 139 162 185 208 231 254 277 300 323

70 93 116 139 162 185 208 231 254 277 300 323
nanos trains /ye ar

Location 11 - SURF DEF

Location 14 - SURF DEF

0.8
0.7

0.8
0.7
0.6

Weighted Frequency

1
0.9

0.6
0.5
0.4
0.3
0.2
0.1
0
24 47

24 47

nanos trains /ye ar

1
0.9

1

70 93 116 139 162 185 208 231 254 277 300 323
nanos trains /ye ar

1
0.9

24 47

24 47

nanos trains /ye ar

1
0.9

1

70 93 116 139 162 185 208 231 254 277 300 323
nanos trains /ye ar

1
0.9

24 47

24 47

nanos trains /ye ar

1
0.9

1

Weighted Frequency

24 47

Weighted Frequency

Weighted Frequency

1

Weighted Frequency

Location 8 - SURF DEF

Weighted Frequency

Weighted Frequency

Location 4 - SURF DEF

0.5
0.4
0.3
0.2
0.1
0
1

70 93 116 139 162 185 208 231 254 277 300 323

24 47

70 93 116 139 162 185 208 231 254 277 300 323
nanos trains /ye arnanos trains /ye ar

Figure 3.10 (continued from previous page): Strain histograms at 14 example locations from surface
deformation.

23

3.2

Strain Rates from Global Positioning System (GPS) Data

Since the early 1990’s, our understanding of plate boundary zone crustal deformation has
been revolutionized by advances in Global Positioning System (GPS) techniques. Amongother things, GPS techniques allow us to directly track the movement of tectonic plates in realtime, measure the rates of crustal deformation within fault systems at plate boundary zones,
and measure the displacement of the Earth’s surface during earthquakes. GPS
measurements are taken at permanent survey points either by intermittent (campaign-style)or continuous collection of phase and pseudorange data from the GPS constellation ofsatellites that orbit the Earth. Using GPS measurements spanning some period of time(usually > 1 year) we can accurately measure (to within a few mm uncertainty) the movementof a survey point on the Earth’s surface relative to another GPS survey point, or relative tosome frame of reference. Such measurements have allowed scientists to determine where
and at what rate tectonic strain is currently accumulating within plate boundary zones (e.g.,
Sagiya et al., 2000; McClusky et al., 2000; Kreemer et al., 2000; Beavan and Haines, 2001).
Although GPS datasets provide a direct measure of contemporary surface strain rates, andcontain important information regarding “rock deformation” hazards, they are only slowly
being introduced in nuclear facility or HLW repository hazard assessments. This isundoubtedly due to the relative newness of the technique.
In most cases, we cannot assume that all of the strain measured by contemporary GPS will
result in an earthquake in the exact spot where the strain is currently accumulating. That is
because some of the strain being measured by GPS is transient, elastic strain that occurs in
the time between large earthquakes on major faults (the interseismic period); this elasticstrain can occur over a wide region on either side of the fault (Figure 3.11a). Moreover, large
displacements of GPS sites during earthquakes (coseismic) or for some period of time
following an earthquake (postseismic) can complicate the interpretation of GPS time series.
Although GPS datasets commonly record transient effects, not representative of the long-term
tectonic deformation, a variety of techniques has been devised to extract the long-term,
tectonic component from the contemporary GPS strain signal. We have developed a method
to extract the long-term tectonic strain signal (relevant to HLW repository siting issues) from a
regional GPS velocity field in the Tohoku region. We also discuss how this tectonic strain
signal might be used in the assessment of future rock deformation/tectonic hazard.

3.2.1

What exactly is GPS measuring?

In the most general sense, when we conduct repeated GPS measurements at a network ofsites, we are measuring the relative movement between GPS survey points; suchmeasurements can then be used to calculate the crustal strain within the GPS network. We
can also estimate the change in position (with time) of a given GPS survey point relative tosome frame of reference; tectonic plates are commonly used to define a reference frame. A
linear fit to timeseries of these positions yields the velocity of the GPS survey point relative
the reference frame.
The motion of tectonic plates relative to some reference frame can be described by the
angular velocity of the plate on a sphere (e.g., Cox and Hart, 1986). The point at which theplate’s angular velocity vector intersects the Earth’s surface is known as a “pole of rotation”.
GPS velocities are often used to estimate the pole of rotation and rotation rates of tectonic
plates (e.g., Beavan et al., 2002; Sella, 2002; Prawirodirdjo and Bock, 2004). Seafloor
spreading rates, oceanic transform fault azimuths and earthquake slip vectors can also be
used to determine tectonic plate motions (e.g., DeMets et al., 1994). Angular velocities oftectonic plates derived from geological data (DeMets et al., 1994) and those derived fromGPS data (e.g., Beavan et al., 2002; Sella et al., 2002; Prawirodirdjo and Bock, 2004) agree
very well in most cases, suggesting that tectonic plate motions have been largely steady overthe last 1-2 Ma, and that velocities of GPS sites within the interiors of large tectonic plates canbe explained by these steady tectonic plate movements. Accurate estimates of tectonic platemovements from GPS can also give us the “plate motion budget” that must be accounted forwithin plate boundary zones when conducting tectonic hazard assessments. The plate motion
budget is the total rate of deformation that must occur across a plate boundary zone, in orderto accommodate the relative motion between the adjacent tectonic plates. We emphasise that24

a key component of any rock deformation hazard assessment is to ensure that all of the platemotion budget in a plate boundary zone is accounted for; in other words, that we know where
all of the rock deformation accommodating the total relative plate motion is occurring.

Figure 3.11: a) An idealized schematic of crustal deformation throughout a seismic cycle, assuming
deformation of an elastic medium. The two blocks are sliding past each other along a strike-slip fault,
and the dashed lines crossing the fault are drawn to illustrate the deformation occurring across the fault
as time progresses. A GPS site is shown for schematic purposes on the left-hand fault block (filled black
ellipse with two lines at the base). (b) Hypothetical displacement of GPS survey point (in a similar
position to that shown in [a]) throughout an earthquake cycle.

GPS velocity and/or strain measurements may also comprise a variety of signals that must be
understood prior to using the GPS measurements for tectonic/rock deformation hazard
assessment, particularly in plate boundary zones. Much of the strain that contemporary GPS
measurements detect is elastic strain that accumulates in the Earth’s crust during the time
between large earthquakes (e.g., McCaffrey et al., 2000; Mazzotti et al., 2000; Wallace et al.,2004). This strain arises when two pieces of crust are prevented from sliding past one anotheralong a fault due to friction on the fault’s surface. The two pieces of crust become stuck
together along the fault, causing the surrounding crust to deform elastically as the crustal
blocks continue to move in the far-field (Figure 3.11a). Throughout this report, we refer to thisphenomenon as interseismic coupling, and the degree of coupling is a quantity that we term
the “coupling coefficient” (φ), which can have a value between 0 and 1. Where φ = 0, the fault
creeps aseismically with no elastic strain accumulation, and for φ = 1, the fault is fully coupled.
Eventually, the stress build-up due to interseismic coupling becomes so large that an
earthquake occurs as the fault fails and the blocks slip rapidly past each other along the fault.
In an earthquake, the elastic strain that has built-up in the surrounding crust during thepreceding interseismic period becomes recovered and is converted to slip across the fault(Fig. 3.11a). In many cases, this elastic strain signal is so large that it can mask strain buildup related to other active faults in the region, particularly at subduction zone margins whereelastic strain due to interseismic coupling on the subduction interface can occur hundreds ofkilometres from the trench (e.g., McCaffrey et al., 2000; Mazzotti et al., 2000; Wallace et al.,2004).

25

Transient deformation related to an earthquake can continue for some period of time following
an earthquake, largely due to: (1) continued slip on the fault or nearby faults (often aseismic,
usually called “postseismic slip”, or “afterslip”) that can occur for a period of a few months toseveral years following the earthquake (e.g., Savage and Svarc, 1997; Heki et al., 1997;
Segall et al., 2000) and (2) viscoelastic relaxation of the mantle and lower crust as it readjuststo the change in the Earth’s crust that has occurred during the earthquake (this can last forseveral decades following the earthquake) (Pollitz, 1997, and references therein). To estimate
the impact of viscoelastic relaxation on the GPS dataset, modelling of this process is often
conducted (e.g., Pollitz, 1997; Pollitz et al., 2006). To correct for coseismic and postseismic
offsets in the GPS timeseries (Figure 3.11b), many studies use timeseries analysis
techniques to estimate the coseismic offset as a step function, while postseismic deformationis commonly approximated by an exponential decay function (e.g., Heki et al., 1997; Hsu etal., 2002; Prawirodirdjo and Bock, 2004). Accounting for and removing the effect of coseismic
and postseismic displacements that occur within a GPS network is critical because these
displacements often produce large, measurable strains, which are much higher than the
yearly-averaged background tectonic strain. Slow slip events (accelerated fault slip eventsoccurring in a matter of days to years, e.g., Hirose et al., 1999, Dragert et al., 2001; Larson etal., 2004) that have been recently documented at several subduction margins must also be
accounted for when assessing the GPS time-series. Transient deformation related to nearby
volcanic activity (e.g., Dixon et al., 1997; Mattioli et al., 1998; Miura et al., 2000) can further
complicate the tectonic interpretation of GPS data, although the influence of volcanic
deformation usually only affects GPS velocities within several km of the volcanic activity.
Seasonal effects are observed in most GPS time-series, and can be estimated and removed(e.g., Nikolaidis, 2002; Williams et al., 2004).

3.2.2

Extracting the component of the GPS signal relevant to rockdeformation hazards

Although numerous types of transient deformation can influence contemporary GPS
measurements, in many cases the influence of these effects on the GPS dataset can be
determined by various modelling techniques. One of the most challenging transient signals todeal with is elastic deformation due to interseismic coupling on major, nearby faults. Strain
rates measured by GPS in zones of active faulting are primarily due to the elastic strainscaused by interseismic coupling on the major faults (Figure 3.11a). Where these elasticstrains occur is not always representative of where it will be converted to permanent rockdeformation (strain) across a discrete fault, and subtle straining due to minor faults (ordistributed strain) may not be easily visible in the presence of the larger elastic strain rates.
To model these elastic strains, we often approximate the Earth’s crust as an elastic halfspace. Measurements of deformation of the crust between earthquakes and during
earthquakes have shown that to first-order, the crust can be approximated as an elastic
material. Mathematical expressions have been derived using elastic dislocation theory to
predict horizontal and vertical ground deformation due to fault slip (dislocations) in an elastic,
half-space (e.g., Mansinha and Smylie, 1971; Okada, 1985, among others).
Some studies suggest that the elastic half-space assumption is too simplistic, and thatlayered space elastic models (e.g., Wald and Graves, 2001; Zhu and Rivera, 2002), ormodels of the crust that have an elastic layer overlying a viscoelastic layer may be more
appropriate (e.g., Thatcher and Rundle, 1984). In particular, elastic, half-space models
perform poorly in comparison to layered space models when dealing with coseismic
displacements due to shallow earthquake ruptures, and in basin environments (the half-spacemodels tend to under-predict the surface displacements in these situations; Wald and Graves,2001). Viscoelastic models are often used because they may help explain discrepancies
between geodetic and geological estimates of fault slip rates that may be due to possible
time-dependent (earthquake cycle-related) deformation (Thatcher and Rundle, 1984; Dixon etal., 2003). Despite some of the obvious limitations of the elastic, half-space approach, it works
remarkably well in characterizing coseismic and interseismic deformation fields measured bya variety of geodetic techniques (e.g., Beanland et al., 1990; Sagiya and Thatcher, 1999;
Mazzotti et al., 2000; Miura et al., 2004b; Nishimura et al., 2004; Meade and Hager, 2005;
McCaffrey, 2005; Wallace et al., 2004, 2007), and is widely-used in crustal deformation
research. For example, the interseismic coupling distribution in the source region of the 2003

26

Tokachi-Oki earthquake and the coseismic slip distribution (both estimated using elastic
dislocation modelling techniques) match well, indicating that the coseismic slip occurred on
the portion of the subduction interface that was previously coupled in the interseismic period
as estimated using elastic, half-space models (e.g., Miura et al., 2004b). Moreover, thecoseismic slip distributions from GPS displacements (using elastic half-space models) inseveral earthquakes coincide well with the region of coseismic slip determined using
seismological data (e.g., Miura et al., 2004b; Subaraya et al., 2006, among others).
Interpretation of GPS velocity fields in plate boundary zones that assume deformation of anelastic medium usually yield an excellent fit to geologically estimated fault slip rates (e.g.,
McCaffrey, 2005; Reilinger et al., 2006; Wallace et al., 2007 among others), with a relatively
small number of exceptions (e.g., Dixon et al., 2003; Bennett et al., 2004; Wallace et al.,2004).
Recently, many studies have shown that GPS velocities measured in zones of active faulting
during the interseismic period are explained by interseismic elastic strains (Figure 3.11a), as
well as long-term rotation of crustal blocks in the deforming zone (e.g., McCaffrey et al., 2000;
McCaffrey, 2002, 2005; McClusky et al., 2001; Meade and Hager, 2005; Wallace et al., 2004,
2007). Methods have been devised by McCaffrey (1995, 2002) and Meade and Hager (2005)
to invert GPS velocities for long-term rotations of tectonic blocks, and elastic strain due tocoupling on block-bounding faults. For the purposes of rock deformation/tectonic hazards
assessment, if the elastic deformation from the inversion (predicted from the fault coupling
parameters) is subtracted from the original GPS velocities, the resulting velocity field will beapproximately free of the elastic effects of interseismic coupling on known, major active faults
in the region. These residual velocity and strain-rate fields may be interpreted in terms ofdeformation due to other minor faults or perhaps zones of distributed deformation.
McCaffrey’s (1995, 2002) method performs a non-linear inversion to simultaneously estimate
the angular velocities of elastic blocks and coupling coefficients on block-bounding faults, togive the best fit to the GPS velocities, and optionally, earthquake slip vectors, and geologicalfault slip rates and azimuths. The data misfit, defined by the reduced chi-squared statistic (χn2),
is minimized. The method also allows us to optimally rotate multiple GPS velocity solutions
into a common reference frame. McCaffrey’s approach also has the benefit of including all of
the bounding plates so that we can establish the plate motion budget that needs to beaccounted for in the plate boundary zone. Once the elastic deformation effects (estimatedusing the elastic block method) have been removed from the GPS velocity field, the residualstrain can be mapped using a variety of methods (e.g., Haines and Holt, 1993; Beavan andHaines, 2001; Sagiya et al., 2000; Miura et al., 2004a). However, we caution that the strainmapping results from GPS should only be used as a rough guide to identify regions where
strain rates are too high for safety. Given the uncertainties inherent in GPS measurements,
GPS techniques are unable to reliably detect strain rates below a certain threshold (thisthreshold will be dependent on the quality of the GPS network and data, and distribution of
GPS sites); thus, GPS should not be used on its own to determine if a location is tectonicallystable (i.e., additional geological and seismological investigations will be needed to confirm
tectonic stability.)

3.2.3

Application of our proposed method to the Tohoku case study

The Tohoku region of Japan provides an excellent location to test our proposed method forusing GPS velocities to determine where tectonic strain will be accommodated by permanent
deformation of the Earth’s crust. Here, we can test whether or not it is possible to isolate
upper plate strain in the face of very large elastic strains due to interseismic coupling on an
offshore subduction margin. In the Tohoku region, subduction of the Pacific Plate beneath
northern Honshu occurs at the Japan Trench (Figure 3.12).
The subduction interface has produced major earthquakes in the past (e.g., Yamanaka andKikuchi, 2004 and references therein), and large portions of the subduction interface are
currently undergoing interseismic coupling. Elastic strains due to this interplate coupling havebeen detected by the continuous GPS network in Tohoku (e.g., Mazzotti et al., 2000;
Nishimura et al., 2004; Suwa et al., 2006). GPS velocities are also influenced by an active
convergent zone in the Sea of Japan, which accommodates underthrusting of the Sea ofJapan beneath the west coast of northern Honshu; some of this convergence may also be27accommodated by faulting onland, near the west coast of northern Honshu (Figure 3.12).
Active reverse faulting occurs within the Backbone Range of the Tohoku region, although thecontraction related to this zone (as measured by contemporary GPS) is largely masked by thecontractional strain from interseismic coupling on the subduction interface (Miura et al.,2004a). Elastic strain due to interseismic coupling on the subduction interface may penetrate
up to a few hundred kilometres inland from the Japan Trench. To isolate the upper plate
tectonic strains relevant to a rock deformation hazard assessment, we apply McCaffrey’s(1995, 2002) elastic block method to the Tohoku region of northern Japan to estimate theelastic strain due to coupling on the Japan Trench and possible elastic effects frominterseismic coupling on offshore faults in the Sea of Japan convergent zone. We thenremove the elastic effects due to interseismic coupling on the subduction zone and theoffshore Sea of Japan convergent zone from the GPS velocities, and map the residual strain.
The assumption is that this residual strain will eventually result in permanent rock deformationin the upper plate, which could pose a problem if a nuclear facility is located in one of thesehigh-strain regions. Our suggested approach of removing ‘known’ elastic strain signals may
also provide a means to facilitate the incorporation of geodetic data into seismic hazard
models. In the following sections, we show an example of the use of this method for Northern
Honshu, while later in the report we discuss the integration of this method into a probabilistic
approach for estimating the rock deformation hazard.

Figure 3.12: Tectonic setting and GPS velocity field in Northern Honshu (shown relative to a fixed
Eurasian plate). Thinner black lines are active fault traces.

28

3.2.4

GPS velocities, uncertainties and logic tree construction

Velocities of GPS sites throughout Japan are derived from combination of SINEX (Solution
INdependent EXchange format) files provided to us by the Geographical Survey Institute(GSI; http://mekira.gsi.go.jp/) from their daily processing of the Geonet Continuous GPS(CGPS) network in Japan (~1200 CGPS sites in total). The daily network processing of theGeonet network is conducted by GSI, with Bernese GPS processing software (Rothacher andMervart, 1996; Beutler et al., 2001) using standard processing methods. The SINEX files we
use are from one day every three months for the period 1996-2004. In order to combine thedaily SINEX files to estimate velocities for the CGPS sites in Japan relative to a known
terrestrial reference frame, we use GLOBK software (e.g., Herring, 2001). To help place theJapanese dataset in a global context, we also use daily solutions from Scripps Institute of
Oceanography processing of the global IGS network of GPS sites (http://sopac.ucsd.edu), as
well as SINEX files from processing of a subset of ~10 Japanese sites and several global
sites that have been submitted by GSI to the Crustal Dynamics Data Information System(CDDIS; http://cddis.gsfc.nasa.gov/). Using GLOBK we estimate a rotation and translation ofeach dataset into the ITRF2000 reference frame (Altamimi et al., 2002), for each day. To
accomplish this, we tightly constrain the coordinates of a subset of the most reliable IGS GPS
stations to their known ITRF2000 values. We do this for each set of daily solutions, obtaining
a time series of site positions in the ITRF2000 reference frame. The ITRF2000 velocities at
each GPS station are calculated by a linear fit to the daily ITRF2000 coordinates. The
uncertainties in the linear fits are derived using a white-noise model, so the uncertainties areseriously underestimated (e.g., Zhang et al., 1997; Williams et al., 2004). We multiply theformal uncertainties by 5 to give “reasonable” values of about 1 mm/year uncertainty inhorizontal velocities for long-running stations within Japan (T. Nishimura, pers. comm.., 2005).
Ideally, the GPS velocity errors should be assessed more rigorously. This will require
maximum-likelihood analysis of (probably daily, perhaps weekly) time series of GPS positions,to define the appropriate noise model for the data and to calculate a realistic velocity
uncertainty (e.g., Williams et al., 2004; Langbein, 2004).
It is very important to avoid the effects on our velocity estimates from earthquakes and slow
slip events. This is more straight-forward in northern Honshu than in other parts of Japan.
For example, slow slip events have not been observed in the CGPS timeseries in northern
Honshu. Large surface displacements were observed at CGPS sites in northern Honshu due
to the 23 September 2003 Tokachi-Oki earthquake and 25 July 2003 northern Miyagi
earthquake. Given that these earthquakes occurred late in the data time series, we omit thedaily solutions from July 2003 and later in the velocity estimation. We also conducted a visualinspection of the daily position time series for all the sites in the dataset to be sure that there
was no non-linear behaviour recorded by the GPS sites that is not representative of steady
movement during the interseismic period. Iwate volcano in central Tohoku experienced
unrest and caused measurable deformation at nearby GPS sites in 1998 (Miura et al., 2000).
We removed the velocities of affected sites from the dataset.
To interpret the GPS site velocities in a tectonically meaningful way, it is important to place
the velocities into a “plate-fixed” reference frame (in this case, we choose Eurasia-fixed). We
do this at the inversion/tectonic modelling stage of this work by estimating a rotation of theentire dataset that minimizes the velocities at sites known to be on the stable Eurasian Plate.
The GPS velocity field in Northern Honshu is shown in a Eurasia-fixed reference frame in
Figure 3.12. In addition to the GPS velocities described earlier, we use published GPS
velocities in the inversion from Heki et al. (1999), Calais et al., (2003), Sella et al. (2002),
Beavan et al. (2002), and Prawirodirdjo et al. (2004). These datasets help us to further place
the Northern Honshu GPS dataset into a regional plate kinematic context. This is critical as
we must be sure to account for all of the possible relative motion between the various tectonic
plates in the system that could be influencing the GPS measurements in Japan.

29

Figure 3.13: Logic Tree (with weightings shown) constructed by expert panel and used for probabilistic
GPS strain map.

For the GPS logic tree (Figure 3.13), we started out with 3 main models. (1) The panel
suggested that we include the residuals to the Nishimura et al. (2004) coupling model as oneof the main branches; this is a good way to explore the influence of an alternative modelling
approach to the elastic block approach (although both approaches do assume deformation ofan elastic medium and use the backslip technique). (2) The second primary branch considers
a situation where Northern Honshu constitutes a single tectonic block (similar to the example
model shown earlier in this report). (3) The third branch includes an additional block boundary
within the Backbone range. We applied a 20% weighting to the first two models, and 60% to
the third branch.
The next set of options branching from the above block model options (2) and (3) that thepanel recommended we include is the possibility of deeper subduction zone coupling. Thisconsists of two branches, one allowing deeper coupling (down to 100 km), with the otheroption limiting the coupling depth to 80 km. We weighted the deep vs. shallower coupling
models equally.
The set of branches following deep/shallow coupling options consider variation in the dip of
the Sea of Japan fault we use, or to run scenarios where coupling on this feature is notincluded at all. We decided on 3 scenarios (1) no interseismic fault coupling in the Sea ofJapan, (2) a Sea of Japan model fault with a 30 degree dip, and (3) the Sea of Japan model
fault has a 50 degree dip. The panel suggested that these should each be weighted equally at
33.3% each.

30

3.2.5

One representative elastic block model for northern Honshu

Model set-up
In Northern Honshu there are two clear tectonic block boundaries, both of which outcrop
offshore: (1) the subduction interface of the Pacific Plate beneath northern Honshu, and (2)the convergent zone in the Sea of Japan. The Japan Trench and Sea of Japan boundaries
are used as the eastern and western boundary of the northern Honshu block, respectively. In
the block modeling, we also define an Okhotsk block (Figure 3.14a), Amurian Plate, Eurasian
Plate, North American Plate, and a Pacific Plate. The existence of these plates has been
proposed by numerous previous studies (e.g., Seno et al., 1996; Wei and Seno, 1998; Heki etal., 1999; Takahashi et al., 1999) and we define the boundaries of these large plates basedon a digital compilation of tectonic plate boundaries by Bird (2003). To define the subduction
interface fault (Japan Trench), we approximate the configurations shown in Hasegawa et al.
(1994) and Mazzotti et al. (2000) (Figure 3.14b). On the subduction interface, we specify
individual nodes defining the interface spaced on average 40 km apart along strike, and at 10km depth intervals between 0 and 110 km (Figure 3.14b). For the purposes of thisdemonstration, we approximate the Sea of Japan convergent zone as a single fault dipping
50° eastward (Figure 3.14b), although alternative dips and fault configurations could be tried
as well. McCaffrey’s (2002) method is used to solve for coupling coefficients at nodes on the
Sea of Japan fault and the Japan Trench.
To represent the change in coupling coefficient (φ) values between adjacent nodes, φ values
on 5x5 km rectangular fault patches between the nodes are estimated by bilinear interpolation.
We also solve for the rotation of the Northern Honshu (NHON), Okhotsk (OKHO), Pacific(PACI), Amurian (AMUR), and North American (NOAM) blocks relative to Eurasia, androtation parameters that rotate each GPS velocity dataset into a Eurasia-fixed reference
frame.
The addition of the larger, surrounding tectonic plates to the model (PACI, AMUR and NOAM)
assists in balancing the plate motion budget. In addition to the GPS velocities we also include
earthquake slip vectors from events on the Japan Trench and on the Sea of Japan
convergent zone (from Harvard CMT, http://www.globalcmt.org/CMTsearch.html) as
these data give us information about the relative motion between the NHON and PACI blocks,
and the NHON and AMUR blocks, respectively. We impose the constraint that the coupling
coefficients at the fault nodes decrease in value down-dip. The down-dip constraint isnecessary when using an elastic dislocation modelling approach for a subduction interface inorder to avoid the unrealistic extensional strain above the up-dip end of the coupled zone thatis predicted by the model when the area up-dip of the coupled zone of a dipping fault slips
aseismically (e.g., McCaffrey, 2002). This means that we cannot discern the actual up-dip
end of coupling on the subduction interface using elastic dislocation methods.
To test if an improvement in fit to the data can be found by not imposing a down-dip decreasein φ, we ran inversions where this constraint is removed. There is no improvement in fit to thedata when we remove this constraint, suggesting that if there are “uncoupled” areas updip of
the “coupled” zone, the data are insensitive to them.
Block model results
The elastic block modelling results are summarized in Figure 3.14. The fit to the GPS
velocities is satisfactory, with χn2 = 2.3, using 1573 data to estimate 230 free parameters. We
estimate significant interseismic coupling on the Japan Trench and in the Sea of Japan,
consistent with previous studies (Mazzotti et al., 2000; Nishimura et al., 2004; Suwa et al.,2006). Our estimates of rotation of the NHON block relative to the PACI and AMUR platessuggests that ~85 mm/year of convergence is accommodated at the Japan Trench, and 15-20mm/year of convergence in the Sea of Japan, yielding a total plate motion budget of ~102-103mm/year that must occur across the PACI/NHON/AMUR plate boundary zone. Some
proportion of this budget could also be accommodated within the Northern Honshu block itself,
and not only on its boundaries, as the elastic block model assumes. The influence on the
GPS velocities from interseismic coupling on the block-bounding faults is shown in Figure31

3.14b. We assume that this component of the GPS velocity field is elastic (temporary)deformation, and will be recovered during slip on the block-bounding faults in futureearthquakes. It is this elastic component of the velocity field (due to interseismic coupling on
the Japan Trench and Sea of Japan faults) that we wish to remove from the GPS velocities inorder to estimate the residual strains which may be due to active faulting within the Tohoku
region.

Figure 3.14: (a) Block boundaries (heavy grey lines) used in the elastic block modelling and fits to the
GPS velocities (red arrows) from the best-fitting model (black arrows). Red arrows on the boundaries
with numbers beside them (mm/year) show the motion of the Northern Honshu block relative to the
Amurian and Pacific Plates. (b) Coupling coefficients on the Japan Trench and Sea of Japan fault
boundaries estimated from the elastic block model. The black vectors show the influence on the GPS
velocity field from the interseismic coupling coefficients on the Japan Trench and Sea of Japan blockbounding faults. Numbers next to the nodes on the Japan Trench and Sea of Japan faults shows thedepth of the fault nodes.

323.2.6

GPS residual upper plate strain map and example location calculations

There are a variety of methods that can be used to convert GPS site velocities to a map ofregional strain. Perhaps the most widely used method is one developed by John Haines andBill Holt (e.g., Haines and Holt, 1993; Beavan and Haines, 2001). To employ their method, agrid (in latitude-longitude) is developed over the area of interest; for the purposes of the
Tohoku case study, we use a ~5x5 km grid. Velocities are modelled as bicubic splines withineach cell, and the inversion attempts both to match the input velocity data and minimize thestrain rates within each cell. A strain-rate variance parameter (1/ν in the terminology ofBeavan and Haines, 2001) is defined, with lower values of 1/ν giving smoother solutions. The
strain-rate variance parameter is chosen such that the sum of squared residuals between themodel and input GPS velocities plus the sum of squares in matching the strain-rate
constraints is approximately equal to the number of degrees of freedom in the GPS data set(i.e., twice the number of velocities). In other words the reduced chi-squared statistic, χ2N =(sum of squared residuals)/(degrees of freedom), is approximately equal to 1.
Figure 3.15 shows the shear and areal strain rates for the GPS velocity field without removing
the elastic strain due to offshore fault coupling. Clearly, these strain rates are quite high (> 80nanostrain/year for areal strain, and > 180 nanostrain/year for the maximum shear strain), but
as discussed throughout this report they are largely elastic strains from interseismic couplingon the major faults (Japan Trench, Sea of Japan) and are unlikely to lead to permanent
deformation in the same location where the strain is currently accumulating. Figures 3.16 and3.17 show the strain field (areal and shear components) after this elastic part of the GPS
velocity field (Figure 3.14b) has been removed from the raw velocity field (Figure 3.12). The
shear strain rate is greatly reduced (Figure 3.17), and there is no coherent pattern of shear
strain (although there are some isolated patches of high shear strain) that might require
further investigation. This is consistent with the lack of evidence for active strike-slip faults
within northern Honshu.
The areal strain is greatly reduced in Figure 3.16 compared to Figure 3.15a. However, azone of elevated contractional strain (~20 nanostrain/year, χ2N = 1.0 case) persists along thecentral Backbone Range, coinciding with a zone of well-documented, active reverse faults
thought to be accommodating a total horizontal slip of ~1-2 mm/year (Active Fault ResearchCenter database; http://www.aist.go.jp/RIODB/activefault/cgi-bin/index.cgi). Our residualstrain rate estimates suggest that ~1-3 mm/year of plate boundary convergence is taken up
within the Backbone Range. Moreover, the Senya Fault (within the Backbone Range; Figure3.12) ruptured in 1896, producing 3.5 m of reverse displacement in some locations (Matsuda
et al., 1980), providing additional evidence for ongoing contraction within the Backbone
Range. Using a backslip approach without block rotations, Miura et al. (2004a) obtained asimilar result of high contractional strain along the Backbone Range when they removed theelastic strain due to subduction coupling (Japan Trench) from the GPS velocity field in
northern Honshu.
In addition to the Backbone Range contraction, we also detect elevated contraction along thewest coast of northern Honshu, possibly due to onshore active faulting associated with the
Sea of Japan convergent zone. However, due to our arguably inaccurate assumption that Sea
of Japan convergence occurs on single dipping fault offshore of the west coast of northern
Honshu (Figure 3.14b) it is possible that we have under-estimated the contraction in the GPS
velocity field due to active faulting onshore near the west coast. For example, much of thecontraction that may be due to active faulting on the west coast has been interpreted aselastic strain due to coupling on the offshore Sea of Japan boundary in our block model. To
incorporate this and other alternative tectonic models into our Tohoku regional strain
estimates, we can implement a probabilistic approach to estimating these strains (seediscussion in next section).

33

Figure 3.15: (a) areal and (b) shear strain rates from the raw GPS velocity field (see Figure 3.14; noelastic strains are removed) estimated using the Haines and Holt method. Strain rates shown innanostrains/year.

34

Figure 3.16: (a) Areal strain in the residual GPS velocity field in Northern Honshu, after the elastic
component of the velocity field from coupling on block-bounding faults (e.g., Figure 3.14b) is removed2
2
from the raw velocity field (χ N = 0.5) (b) Areal strain in the residual GPS velocity field for χ N = 1.0 (c)2Areal strain in the residual GPS velocity field for χ N = 2.0. Strain rates shown in nanostrains/year. Notethe colour scales are identical to those in Figure 3.15a, and the contour intervals are 40 nanostrains.

35

Figure 3.17: (a) Shear strain in the residual GPS velocity field in Northern Honshu, after the elastic
component of the velocity field from coupling on block-bounding faults (e.g., Figure 3.14b) is removed
from the raw velocity field (χ2N = 0.5) (b) Shear strain in the residual GPS velocity field for χ2N = 1.0 (c)2Shear strain in the residual GPS velocity field for χ N = 2.0. Strain rates shown in nanostrains/year.
Note the colour scales are identical to those in Figure 3.15b, and the contour intervals are 40
nanostrains.

In Figure 3.18 we show the strain rate histograms for all of the various GPS models
developed. The plots reveal a very wide range of strain rates at many locations. These ranges
reflect the variety of parameter choices between models. Bi-modal plots often result fromchoices made regarding locking depths of subduction faults at the Japan trench and on thefault system along the west part of Honshu and offshore. As with the surface deformation
dataset, there are marked differences in strain rates between some of the example locations,suggesting the GPS data does reflect regional variation in hazard. At many of the example
locations the strain rates are in the range of 5-50 nanostrains/year, similar to the surface
deformation strains. We discuss comparisons between datasets in Section 3.4.

36

Location 5 - GPS

Location 1 - GPS
0.4

0.4

0.35

0.35
0.3
Weighted Frequency

Weighted Frequency

0.3
0.25
0.2
0.15
0.1
0.05

0.25
0.2
0.15
0.1
0.05

0

0
1

5

9 13 17 21 25 29 33 37 41 45 49 53 57 61 65 69 73 77 81 85 89

15

9 13 17 21 25 29 33 37 41 45 49 53 57 61 65 69 73 77 81 85 89

nanotrains/year

nanotrains/year

Location 6 - GPS

Location 2 - GPS
0.4

0.4

0.35

0.35
0.3
Weighted Frequency

Weighted Frequency

0.3
0.25
0.2
0.15
0.1
0.05

0.25
0.2
0.15
0.1
0.05

0

0
1

5

9 13 17 21 25 29 33 37 41 45 49 53 57 61 65 69 73 77 81 85 89

1

5

9

13 17 21 25 29 33 37 41 45 49 53 57 61 65 69 73 77 81 85 89

nanotrains/year

nanostrains/year

Location 3 - GPS
0.4
0.35

Weighted Frequency

0.3
0.25
0.2
0.15
0.1
0.05

89

85

81

77

73

69

65

61

57

53

49

45

41

37

33

29

25

21

17

9

13

5

1

0
nanotrains/year

Location 4 - GPS

Location 8 - GPS

0.4

0.4

0.35

0.35
0.3

0.2
0.15

89

85

81

77

73

69

65

61

57

53

49

45

41

37

33

29

25

21

17

1

89

85

81

77

73

69

65

61

57

53

49

45

41

37

33

29

25

21

17

13

0

9

0
5

0.05

9

0.1

0.05

13

0.1

0.25

5

Weighted Frequency

0.2
0.15

1

Weighted Frequency

0.3
0.25

nanostrains/year

nanotrains/year

Figure 3.18 (continued on next page): Strain rate histograms derived from GPS at example locations
1-8. See Fig 3.1 for locations, and caption of Figure 3.8 for full explanation of construction of histograms.
Note: there are no GPS data for Location 7. Data for all the locations are presented in the Appendix.

37

Location 12 - GPS

Location 9 - GPS
0.4

0.4

0.35

0.35
0.3
Weighted Frequency

0.2
0.15
0.1
0.05

0.25
0.2
0.15
0.1
0.05

0

Location 10 - GPS

89

85

81

77

73

69

65

61

57

53

49

45

Location 13 - GPS

0.4

0.4

0.35

0.35
0.3
Weighted Frequency

0.3
Weighted Frequency

41

nanostrains/year

nanostrains/year

0.25
0.2
0.15
0.1

0.25
0.2
0.15
0.1
0.05

0.05

0

0
1

5

9

1

13 17 21 25 29 33 37 41 45 49 53 57 61 65 69 73 77 81 85 89

59

13 17 21 25 29 33 37 41 45 49 53 57 61 65 69 73 77 81 85 89
nanostrains/year

nanostrains/year

Location 14 - GPS

Location 11 - GPS
0.4

0.4

0.35

0.35
0.3
Weighted Frequency

0.3
Weighted Frequency

37

33

29

25

21

17

13 17 21 25 29 33 37 41 45 49 53 57 61 65 69 73 77 81 85 89

9

9

1

5

13

0

1

5

Weighted Frequency

0.3
0.25

0.25
0.2
0.15
0.1

0.25
0.2
0.15
0.1
0.05

0.05

0
1

89

85

81

77

73

69

65

61

57

53

49

45

41

37

33

29

25

21

17

9

5

13

1

0

5

9 13 17 21 25 29 33 37 41 45 49 53 57 61 65 69 73 77 81 85 89
nanostrains/ye arnanostrains/year

Figure 3.18: continued from previous page.

3.3

Strain Rates from Seismicity Data

Our methodology to quantify and map seismic strain rates is to first develop a seismicity
model in the same overall way that background seismicity models are developed for PSHAs.
We apply two different methods to develop the seismicity model, each providing a set of
gridded point sources that have a set of earthquake parameters assigned to them, and withthe parameters described by the Gutenberg-Richter relationship. The first is the ‘traditional’
method of defining large seismicity zones from seismotectonic considerations, assigningseismicity parameters to each zone, and then uniformly distributing the parameters across thezone. The second method is to allow the seismicity parameters to vary within the zones
according to the spatial distribution of seismicity within the zone. The latter method wasdeveloped by Frankel (1995) to characterise the PSH from background earthquakes, and has
been successfully adapted to develop PSH models in other regions (e.g. Stirling et al. 2002).
Use of the two alternate methods for treatment of the seismicity data represents important
epistemic uncertainties in the seismic strain model. The traditional method allows for the
possibility that the current seismicity patterns (varying over distances of tens of kilometres) do
not necessarily represent the long-term parent distribution of seismicity, and instead relies on
the broader seismotectonic regions as the guiding definition. The spatially-varying method incontrast assumes that the current seismicity patterns are a reasonable representation of longterm seismicity. Frankel (1995) initially made such arguments based on a spatial correlation
between small, frequently occurring events and less frequently-occurring magnitude 5 and38

larger earthquakes.
These alternative modelling techniques, along with alternative
regionalisation (area source zonation) schemes, use of declustered versus raw seismicity
catalogues, and the uncertainties in the seismicity parameters are the basis for definition oflogic trees for the seismic strain rate model.

3.3.1

Seismicity Dataset

The source of the seismicity data used to develop seismic strain rates in this study is theJapanese Meteorological Agency (JMA) network catalogue. Earthquakes have been routinely
recorded since 1926 at a detection threshold magnitude of around 4.5 (e.g. Stirling et al.
1996). That detection threshold has improved with time and progressive development of the
JMA network (Figure 3.19). Our catalogue is also supplemented by a complete record ofmajor historical earthquakes (M>6.9) for the period 1581 to 1925 (Wesnousky et al. 1984)(Figure 3.19).

Figure 3.19: Major historical and recent earthquakes (large black circles and dates showing M>7
epicentres - Wesnousky et al., 1982) and smaller magnitude instrumental seismicity (grey small symbols- JMA catalogue from 1926-2006) in Tohoku

These alternative modelling techniques, along with alternative regionalisation (area source
zonation) schemes, use of declustered versus raw seismicity catalogues, and theuncertainties in the seismicity parameters are the basis for definition of logic trees for theseismic strain rate model.

3.3.2

Logic Tree Structure

Component parts of the logic tree (Figure 3.20) include assessment of the seismicity
parameters (a-value and b-value of the Gutenberg-Richter relationship; logN=a-bM ),calculated by the maximum-likelihood method of Weichert (1980), which allows the use ofdifferent magnitude completeness levels for various time periods to calculate parameter b.
The Mmax is based on the maximum magnitudes assigned to the various parts of the Tohoku

39

region by the Research Group for Active Faults of Japan (1991), which is based on a
combination of historical earthquake information and consideration of the likely maximum
magnitudes derived from the fault database. The final step in the seismic strain model is toconvert the seismicity rates into equivalent strain rates. This is achieved by the method ofKostrov (1974), which converts the equivalent seismic moment rate from each of the
Gutenberg-Richter-distributed earthquakes into strain rate through assumption of a crustal
volume. Our convention is to use 5x5 km area and 20 km thick, (20 km being the average
depth to the base of crustal seismicity), based on interpretation of seismicity cross sections.

Catalogue

Zonation

Modelling
Method0.3

0.4

3 zonesContinuously
Varying
Seismicity

a-valueb-valuea-v alue 0.2
30km gaussiansmoot hing k ernalsigma
a-v alue

0.6
15km gaussiansmoot hing k ernal
sigma

Seismicity
rates tostrainrates via
Kostrov
equationa-v alue
0.8

Raw

0.2
10km gaussiansmoot hing k ernal
sigma

0.2

Mmax

1 large zone
0.4

Uniform
Area source

4 zones
0.7
Sigma =
sqrt(N)

0.2

DeclusteredSigma =
b/sqrt(N)Sigma =
0.2 for all zones
Except0.3 for 1 large
zone

Figure 3.20: Logic tree structure developed for the seismicity strain model. This incorporated the expert
elicitation meeting results of July 2005. The decimals shown are the weights assigned to each branch,and the schematic bell curve indicates that the parameter has a continuous normal distribution.

The logic tree considers both raw and declustered (i.e. aftershocks removed) catalogues,
because preliminary seismicity strain models produced strain rates an order of magnitude
lower than the rates derived from the geodetic model, leading to the conclusion that the“Gardner and Knopoff” declustering programme was too severe in its identification andsubsequent removal of aftershocks. Three source zonation schemes were developed (Figure3.21).
Weights are given to all of the branches of the logic tree, incorporating expert opinion andconfidence in the use of a raw versus declustered catalogue, choice of zonation scheme, theuse of continuously-varying versus uniform area sources approaches for parameterisation ofseismicity, the range of seismicity parameters (the a and b-values of the Gutenberg-Richter
relationship) and the choice of maximum magnitude for the respective regions. The seismicity
parameter “branches” were represented by continuous distributions, whereby the choice ofvalues of a given parameter are an infinite number within minimum and maximum bounds,
rather than a few discrete choices (Figure 3.20).

40

Figure 3.21: The seismicity of Tohoku (JMA catalogue, 1926 – 2004, for M>3 and 0 to 20 km depth),
and the three alternative seismic zonation schemes used in the seismicity strain rate model. The
seismicity parameters are calculated by way of the maximum-likelihood method of Weichert (1980).
Zonation scheme 1 comprises the following: Sea of Japan, West Honshu, and All of East Honshu (i.e.
combining Onshore and Offshore East Honshu zones). Zonation scheme 2 comprises: Sea of Japan,
West Honshu, Onshore East Honshu, and Offshore East Honshu, and Zonation scheme 3 comprises:
All zones combined into one large zone.

3.3.3

Results and Discussion

A seismicity-derived strain rate map, comprising the highest weighted branches of the logic
tree is shown in Figure 3.22. It shows considerably higher rates of activity along the JapanSea coast (1993 major aftershock sequence) than the rest of the region, and somewhat
higher seismicity in the Ou-Backbone Range region, and in a small part of the southern
Kitakami region. This last area is somewhat anomalous compared with other strain models,and may be the result of the source area zonation scheme in the particular model. In Figure
3.23 we show histograms of strain rates for the 14 example location derived by Monte Carlosampling of a large number of branches of the logic tree. For each of example locations
distributed throughout Tohoku (including a Sea of Japan location), 1000 samples ofparameters were taken from the logic tree for each of the locations, and 1000 strain rates
calculated from each set of parameter values. The 1000 strain rates were then combined to
produce the histograms shown in Figure 3.23. The typical range of strain rates observed at
the 14 example locations is 1-50 nanostrains/year with a mode toward the lower end (Figure3.23).
Further consideration of the completeness of the seismicity catalogue is in order. There is atendency for seismicity to cluster in regions of volcanism, but at low to moderate magnitudelevels, with rare large earthquakes, usually located between the volcanic clusters. We haveexamined focal mechanisms of earthquakes to see if normal mechanisms occur nearvolcanoes consistent with extension due to magma intrusion, but have found no correlation.
Essentially all focal mechanisms are reverse or strike-slip. Our analysis supports the41

conceptual model of Hasegawa et al. (2000) (Figure 3.24), which suggests that small tomoderate magnitude earthquakes are concentrated near the thermally weakened crust of the
volcanic clusters, and large magnitude earthquakes occur rarely in the relatively stronger
crust between volcanic clusters. A further iteration of the seismicity methodology would be to
distinguish between tectonically and volcanically derived seismicity in the source zonation andto examine seismicity parameters on this basis.

Figure 3.22: Seismicity strain map constructed from c. 400 year long record of historical andinstrumental seismicity according to the logic tree shown in Figure 3.20. This plot is for the highest
weighted realisation of the logic tree.

42

Location 5 - SEIS

0.4

0.4

0.35

0.35
Weighted Frequency

Weighted Frequency

Location 1 - SEIS

0.3
0.25
0.2
0.15
0.1
0.05
0

0.3
0.25
0.2
0.15
0.1
0.05

1
5
9
13
17
21
25
29
33
37
41
45
49
53
57
61
65
69
73
77
81
85
89

1
5
9
13
17
21
25
29
33
37
41
45
49
53
57
61
65
69
73
77
81
85
89

0

nanos trains /ye arnanos trains /ye ar

Location 6 - SEIS

0.4

0.4

0.35

0.35

0.3

Weighted Frequency

Weighted frequency

Location 2 - SEIS

0.25
0.2
0.15
0.1
0.05
0

0.3
0.25
0.2
0.15
0.1
0.05
1
5
9
13
17
21
25
29
33
37
41
45
49
53
57
61
65
69
73
77
81
85
89

1
5
9
13
17
21
25
29
33
37
41
45
49
53
57
61
65
69
73
77
81
85
89

0

nanos trains /ye arnanos trains /ye ar

Location 7 - SEIS

0.4

0.4

0.35

0.35
Weighted Frequency

Weighted Frequency

Location 3 - SEIS

0.3
0.25
0.2
0.15
0.1
0.05
0

0.3
0.25
0.2
0.15
0.1
0.05

1
5
9
13
17
21
25
29
33
37
41
45
49
53
57
61
65
69
73
77
81
85
89

1
5
9
13
17
21
25
29
33
37
41
45
49
53
57
61
65
69
73
77
81
85
89

0

nanos trains /ye arnanos trains /ye ar

Location 8 - SEIS

0.4

0.4

0.35

0.35
Weighted Frequency

Weighted Frequency

Location 4 - SEIS

0.3
0.25
0.2
0.15
0.1
0.05
0

0.3
0.25
0.2
0.15
0.1
0.05

1
5
9
13
17
21
25
29
33
37
41
45
49
53
57
61
65
69
73
77
81
85
89

1
5
9
13
17
21
25
29
33
37
41
45
49
53
57
61
65
69
73
77
81
85
89

0

nanos trains /ye arnanos trains /ye ar

Figure 3.23 (continued on next page): Strain rate histograms derived from seismicity at 14 example
locations. See Figure 3.1 for locations, and caption of Figure 3.8 for full explanation of construction ofhistograms.

43

Location 12 - SEIS

0.4

0.4

0.35

0.35

0.3

Weighted Frequency

Weighted frequency

Location 9 - SEIS

0.25
0.2
0.15
0.1
0.05
0

0.3
0.25
0.2
0.15
0.1
0.05

1
5
9
13
17
21
25
29
33
37
41
45
49
53
57
61
65
69
73
77
81
85
89

1
5
9
13
17
21
25
29
33
37
41
45
49
53
57
61
65
69
73
77
81
85
89

0

nanos trains /ye arnanos trains /ye ar

Location 13 - SEIS

0.4

0.4

0.35

0.35
Weighted Frequency

Weighted frequency

Location 10 - SEIS

0.3
0.25
0.2
0.15
0.1
0.05
0

0.3
0.25
0.2
0.15
0.1
0.05

1
5
9
13
17
21
25
29
33
37
41
45
49
53
57
61
65
69
73
77
81
85
89

1
5
9
13
17
21
25
29
33
37
41
45
49
53
57
61
65
69
73
77
81
85
89

0

nanos trains /ye arnanos trains /ye ar

Location 14 - SEIS

0.4

0.4

0.35

0.35
Weighted Frequency

Weighted Frequency

Location 11 - SEIS

0.3
0.25
0.2
0.15
0.1
0.05
0

0.3
0.25
0.2
0.15
0.1
0.05

1
5
9
13
17
21
25
29
33
37
41
45
49
53
57
61
65
69
73
77
81
85
89

1
5
9
13
17
21
25
29
33
37
41
45
49
53
57
61
65
69
73
77
81
85
89

0

nanos trains /ye arnanos trains /ye ar

Figure 3.23 (continued from previous page)

Figure 3.24: Schematic illustration of across-arc vertical cross-section of the crust and uppermost
mantle, showing the characteristic shallow seismic activity and deformation of the crust beneath NE
Japan (From Hasegawa, 2000).

44

3.4

Strain Rate Comparisons

Much of the work of the Tohoku case study has thus far gone into separately developing thethree strain rate models and developing a methodology. Further work is required to develop
the best approaches to integrating the three measures of strain, given the various strengths
and weaknesses of the three datasets (see earlier sections), along with the underlying
assumptions and uncertainties regarding what the individual datasets represent. The GPS
dataset reflects the spatial and temporal accumulation of elastic strain in the crust, whereas
the seismicity and geological data reflect the release of elastic strain across the region overtime. In Figures 3.6, 3.7, 3.8, 3.9, 3.10, 3.16, 3.18, 3.22 and 3.23, strain maps for the whole
region and strain rate estimates for example locations are fairly similar.

Surface Deformation
GPS Areal

Figure 3.25: Comparison of “bestestimate” strain rate maps derived fromeach of the three rock deformation
datasets. Note that the surface
deformation map covers only the
Tohoku area. Also, note that the GPS
strain map represents the weighted
average of all the strain maps in the
GPS logic tree, not just the highest
weighted strain map. To allow easiercomparison with the seismic strain map,the GPS map takes the absolute valueof strain, with no distinction between
extension and contraction.

Seismicity

45

A reasonable correlation is evident from the seismic, surface deformation, and GPS strain
rate histograms and maps, particularly in terms of the roughly north-south trending zone ofhigh strain rate along the central mountains of Tohoku. Higher strain rate zones are also
evident in all three datasets in the west, and there are lower strain rates in the east of thestudy area (Figure 3.25). Therefore, it appears that the strain accumulation measured by GPSis largely elastic and is representative of strain that will be released in future by tectonic
processes.
This “first-order” result is very encouraging for using the rock deformation datasets in acomprehensive manner to estimate future rock deformation hazards. The general similarity inresults between strain accumulation and strain release signatures also means that we can
combine strain rate signatures by way of a weighted averaging procedure to provide a single
measure for comparison of example locations. We should take care however to explore thedifferences between strain models to understand whether differences have a physical basis or
are a consequence of modelling decisions.
When comparing the strain rates seeking similarities and differences, it is also very important
to consider uncertainties is the strain rate values. The probabilistic method allows for formal
assessment of uncertainty in the strain rate estimates. Uncertainty analysis has not yet been
fully explored but in Figure 3.26 we show mean and one standard deviation ranges for strain
rate at four example locations. This figure illustrates general agreement in the strain rates
from different methods (the very high surface deformation rate at location 8 is the result ofvery close proximity to an active fault), and consistent differences in strain rate and thus likely
future hazard between some of these example locations.

Figure 3.26: Strain rate uncertainties at four example locations. The weighted average rate (WGT) is
derived from an equally weighted average of each of the GPS, surface deformation (SURF DEF), andseismicity (SEIS) values.

46

3.5

From Strain Rate to Hazard

Seismological and geological data are important inputs into probabilistic seismic hazard
assessments (e.g., Reiter, 1990; Senior Seismic Hazard Analysis Committee, 1997; Frankel
et al., 2002; Stirling et al., 2002), which have been widely used to estimate the seismic
shaking hazard at many locations around the world, including detailed PSHA work on somekey nuclear facilities (e.g., Stepp et al., 2001).
Methodologies have been developed for PSHA that can embrace a variety of different
techniques (and/or tectonic models) that might be used to convert seismological andgeological data into a probability of “seismic hazard”.
Hazard estimates using a variety of techniques and datasets are combined, by weighting theresults of different models, data, and/or expert opinion in a logic tree environment, which isthen sampled (thousands of times) by Monte Carlo sampling methods to provide a probability
distribution of resulting seismic hazard estimates that collectively expresses the epistemic
uncertainty.
Here, we use a probabilistic method to assess the “rock deformation hazard” (i.e., likelihoodof permanent rock deformation occurring at a given point on the Earth’s surface in the future)by using GPS, geological and seismological data to obtain contemporary surface strain rates.

3.5.1

Hazard Calculation

In Sections 3.2, 3.3, and 3.4 we have presented results in terms of discrete distributions ofannual strain rate across the region or at particular example locations. In this section we
explore the concept of quantifying the probability of exceedance of a given annual strain rate,on the basis that it will eventually be possible to specify what strain rate could represent a
hazard to the viability of a HLW repository.
This discussion is, of course, preliminary because repository concepts are not fully developed
at this stage and the nature and magnitude of significantly damaging strains will vary withdesign. It should be noted that, although a repository may sustain damage, this might not leadto unacceptable radionuclide releases or radiological impacts. We demonstrate how we candevelop the discrete distributions of annualised strain rate into the probability of exceedance
for the same suite of strain rates (i.e. cumulative probability), which are analogous to hazard
curves if strain rate is considered a proxy for tectonic hazard. These are produced byprogressively summing up the strain rate probabilities from each of the bins in the discrete
histograms from highest to lowest strain rate.
The full distribution of strain rate uncertainty seen in the discrete histograms is therefore
incorporated into the cumulative probability plots and controls the shape of the curve. The
probability of exceedance for the suite of strain rates is therefore controlled by the distributionof strain rate uncertainty in the discrete histograms.
The cumulative probability histograms (Figure 3.26) show the probability of exceedance forthe same suite of strain rates shown in the discrete histograms. On each of these plots weshow an equally weighted average model, and in Figure 3.28 we show the suite of equallyweighted models for all 14 example locations. This plot indicates (without taking account of
uncertainty) a wide range in probability of exceedance of medium strain rates of say, 10
nanostrain/year.
At some near-fault locations such as 7, 2 and 8 there is a >60% probability of exceeding sucha medium strain rate, while on elevated plateaux associated with volcanic clusters and distantfrom active faults, the probability that the annual strain rates exceeds 10 nanostrains/year is<20%. In the region east of the volcanic arc the probability of exceeding 10 nanostrains/year
is also relatively low, about 40%.

47

Location 5

Location 1

1

1
GPS

SEIS

Probability of Exceedance

Probability of Exceedance

GPS
0.8

SURF DEF
EQ W GTD AV

0.6
0.4
0.2
0
1

10

100

SURF DEF

0.8

SEIS
EQ W GTD AV

0.6
0.4
0.2
0

1000

1

nanostrains/year

Location 2

1

10nanostrains/year

100

Location 6

1

GPS

SURF DEF

Probability of Exceedance

Probability of Exceedance

GPS
0.8

SEIS
EQ W GTD AV

0.6

1000

0.4
0.2
0

SURF DEF

0.8

SEIS
EQ W GTD AV

0.6
0.4
0.2
0

1

10

100

1000

1

nanostrains/year

10

Location 3

nanostrains/year

100

1000

Location 7

1

1
Probability of Exceedance

Probability of Exceedance

SEIS

GPS

0.8

SURF DEF
SEIS

0.6

EQ W GTD AV
0.4
0.2
0

SURF DEF

0.8

EQ W GTD AG
0.6
0.4
0.2
0

1

10

nanostrains/year

100

1000

1

10

100

Location 8

Location 4
1

1

GPS

GPS
SURF DEF

0.8

Probability of Exceedance

Probability of Exceedance

1000

nanostrains/year

SEIS
EQ W GTD AV

0.6
0.4
0.2
0

SURF DEF

0.8

SEIS
EQ W GTD AV

0.6

0.4

0.2

0

1

10

nanostrains/year

100

1000

1

10

100
nanos trains /ye ar

Figure 3.27: (continued on next page): Cumulative strain rate plots for the 14 example locations. EQ
WGTD AV is equal weighted average

48

1000

Location 9

1

GPS

0.8

Probability of ExceedanceGPS
Probability of Exceedance

Location 12

1

SEIS
EQ W GTD AV

0.6
0.4
0.2

0.8

SURF DEF
SEIS

0.6

EQ W GTD AV

0.4
0.2
0

0
1

10

nanostrains/year

100

1

1000

10

Location 10
1

GPS

GPS
Probability of Exceedance

SURF DEF

0.8

SEIS
EQ W GTD AV

0.6
0.4
0.2

SURF DEF

0.8

SEIS
EQ W GTD AV

0.6
0.4
0.2
0

0
1

10

nanostrains/year

100

1

1000

10

100

1000

nanostrains/year

Location 11

1

Location 14
1
GPS

GPS

SURF DEF

0.8

Probability of Exceedance

Probability of Exceedance

1000

Location 13

1

Probability of Exceedance

100
nanostrains/year

SEIS
EQ W GTD AV

0.6
0.4
0.2
0
1

10

100

0.8

SURF DEF
SEIS

0.6

EQ W GTD AV

0.4
0.2
0

1000

1

nanostrains/year

10

nanostrains/year

Figure 3.27: (continued from previous page)

49

100

1000

Figure 3.28: Plot of probability of exceedance of weighted average annual strain rates for 14 example
locations. Example locations are approximately classified into tectonic domains of the study region. See
Figure 3.1 for location of the 14 example locations50

3.5.2

Probability of Exceedance of Threshold Strain Rates

Discussion in this section demonstrates how strain rates could be assessed against theprobability of exceeding defined thresholds. The procedure is demonstrated by evaluating theprobability of exceeding a strain rate equivalent to 2 m of deformation in a 5x5 km site area,
within a range of time periods. Deformation in the rock mass within the 5 x 5 km site area
could be taken up in a number of ways: (a) concentrated along existing known faults (Figure3.29 A), (b) concentrated along fractures and faults that represent initial development of a
fault at the site (Figure 3.29 B); (c) displacement along principal displacement zones and theassociated process zone of an un-recognised active fault (Figure 3.29 C); (d) uniformly across
many small fractures (Figure 3.29 D). For example, two metres of deformation is consistent
with the localised surface fault rupture associated with an earthquake a little larger thanmagnitude 7 in the Tohoku region (Wesnousky et al, 1984; Matsuda, 1977), rather similar to
the 1896 Rikuu earthquake on the Senya Fault. Two metres of displacement might also
represent repeated surface rupture associated with several earthquakes in the M 6.5-7 range.
Depending on how the deformation is distributed and on repository design, this deformation
may or may not pose problems to repository integrity.
Repository design concepts vary in their tolerance to direct tectonic displacement. The
tolerance of slim (e.g. 1 m diameter and several metres long) canisters may be as little as 0.1m, while large caverns with thick backfill may be able to accommodate as much as a metre or
so of displacement. Therefore siting confidence, and the probability of exceeding a thresholdstrain value that constitutes possible difficulties with the site, will vary with repository concept.
For example, if the repository has a high level of tolerance to deformation then, underscenario (d) above, the tolerable overall strain rate may be far in excess of the 2 m example
value we use here in this case study. Under scenario (b) and (c) newly establishing fault
zones often have a wider process zone than well-established faults (Tchalenko andAmbraseys, 1970; Sylvester, 1988; Boullier et al., 2004; Katz et al, 2004), so more of the 5x5
km site might be subject to deformation. If, under scenario (a), deformation is concentrated ina narrow zone, then much of the 5x5 km site area would not be subject to the deformation.
Sibson (2003; and references therein) has examined the worldwide data for the width of
deformation in historic surface rupturing earthquakes. Variations depending on fault style,geometric complexity, and fault evolution have been noted. This relationship applies to single
events at a site and does not take into account possible fault evolution over periods such as 1
million years when recurrence rupture at the same location is likely.
For this Case Study, we evaluate the probability that 2 m of fault displacement could occur
within a 5x5 km site area within periods of 10,000 years, 100,000 years, and 1 million years.
Two metres of displacement across a 5 km wide zone in 10,000 years equates to a strain rateof 40 nanostrains/year. In 100,000 years the strain rate must exceed 4 nanostrains/year, andmust exceed 0.4 nanostrains/year in 1 million years. Based on the discussion above, these
rates may well be quite conservative for most cases. The combination of a repository designthat is not able to sustain very much deformation coupled with a wide deformation zone
associated with a newly developing fault represents the strain threshold case considered inthis report. We present the analysis here merely to demonstrate the application of the method.
The analysis will need to be iteratively revised as new databases are developed and therepository design concept is further developed.
As noted above, even 2 m of localised strain along a single major discontinuity (fault plus
immediate process zone) that transects a small region of a repository might not lead to
unacceptably increased radionuclide releases. It is an important aspect of probabilistic
methodology that probability is combined with calculated quantitative consequences toestimate hazard. This can be done in the formal sense of a calculated risk or as adisaggregated presentation of consequence and likelihood of occurrence in different time
periods, as is currently being discussed in the development of regulatory standards in Japan.
Our analysis of the region and 14 example locations with respect to the probability ofexceeding nominal (example) threshold values of strain of 40, 4, and 0.4 nanostrains/year,
over periods of 10,000 years, 100,000 years and 1 million years respectively, is illustrated in
Figures 3.30 and 3.31. We emphasise our analysis demonstrates both the general internal51

consistency of our probabilistic methods across the rock deformation datasets and themethod that could be used to make comparative assessment between sites. The analysis
suggests that for our example strain threshold value there is significant variability between
locations for the 10,000 year period and (to a lesser extent) for the 100,000 year period. At 1
million years, the likelihood of exceeding the example value is uniformly close to or equal toone, across the whole region. We emphasise again that these probabilities cannot be used on
their own as criteria to define site suitability. They must be combined with host rock anddesign specific consequence analyses to assess whether such strains would have any
significant impacts.
We expect that, as the databases underpinning the hazard continue to be developed,
increased variability in deformation hazard will emerge. At present we use regional coverages
that tend to smear the strain across broad areas; when analysed at a finer scale, we would
expect there to be more variability in strain rate magnitudes than is evident from the regional
scale maps. Strain rates derived from GPS, seismicity and surface tilting vary only slowly
across the study region, but much of the regional strain will, in time, localise to future
narrower zones of deformation (faults).

Figure 3.29: Schematic illustration of modes whereby strain might occur in or around a repository:
A: all along major faults (old or new) some distance away
B: distal strain uptake from A on small (c.100 m radius) fractures in the repository
C: hidden/new active fault below the repository
D: volume deformation along many small, distributed fractures/planes

52

SURFACE DEFORMATION

Figure 3.30a: Probability ofexceedance of 40
nanostrains/year (10,000
year criterion) for each rockdeformation dataset

53

SURFACE DEFORMATION

Figure 3.30b: Probability ofexceedance of 4 nanostrain/year(100,000 year criterion) for eachrock deformation dataset

54

SURFACE DEFORMATION

Figure 3.30c: Probability ofexceedance of 0.4 nanostrain/yr(1,000,000 year criterion) for eachrock deformation dataset

55

40 nanostrains/yr (~2 m in 10 kyr)
LocationGPS

SURF DEFSEIS

WGTD AV1

6.1E-020

2.2E-03

2.1E-022

8.9E-021

3.7E-03

3.6E-013

1.2E-01

2.0E-03

1.1E-02

4.3E-024

5.7E-020

5.6E-01

1.7E-035

6.1E-020

2.0E-03

2.1E-026

3.5E-0104.4E-03

1.2E-017NA

9.5E-02

9.0E-01

5.0E-018

1.2E-011

1.4E-01

4.2E-019003.7E-03

1.8E-0310

5.3E-020

1.5E-03

1.8E-0211

1.8E-01

1.3E-02

1.4E-02

6.9E-021200

8.3E-03

2.8E-0313

6.1E-02

4.5E-01

1.3E-03

1.7E-0114

2.1E-020

1.6E-03

7.6E-03

4 nanostrains/yr (~2 m in 100 kyr)
LocationGPS

SURF DEFSEIS

WGTD AV1

8.1E-011

5.0E-01

7.7E-012

9.7E-011

5.7E-01

8.5E-013

9.2E-01

8.3E-01

6.4E-01

7.9E-014

6.0E-010

7.6E-01

3.9E-015

7.4E-011

5.1E-01

7.5E-016

9.7E-011

6.2E-01

8.6E-017NA

9.5E-02

9.0E-01

5.0E-018

7.1E-011

8.2E-01

8.4E-019

7.5E-010

6.4E-01

6.9E-0110

8.2E-011

3.6E-01

7.3E-0111

8.9E-011

6.4E-01

8.5E-0112

3.8E-011

8.8E-01

7.5E-0113

8.9E-011

5.1E-01

8.0E-0114

6.9E-011

4.8E-01

7.2E-01

0.4 nanostrains/yr (~2m in 1 Myr)
LocationGPS

SURF DEFSEIS

WGTD AV

1

1

1

1

1

2

1

1

1

1

3

1

1

1

1

4

1

1

1

1

5

1

1

1

1

6

1

1

1

1

7

NA

1

1

1

8

1

1

1

1

9

1

1

1

1

10

1

1

1

1

11

1

1

1

1

12

1

1

1

1

13

1

1

1

1

14

1

1

1

1

56

Figure 3.31: Probability of exceeding
an example specified threshold strain
value equivalent to 2 m of strain
distributed across a site or localised
on a major discontinuity at the 14
example locations.

4 Volcanic Hazard Analysis
4.1

Introduction to the methodology

Worldwide, formal methods are emerging to assess volcanic hazards at sites of a variety offacilities, such as nuclear power plants and geological repositories for radioactive wastes. For
example, the International Atomic Energy Agency is now revising guidelines for theassessment of volcanic hazards at nuclear facilities. Volcanic hazard assessments have been
completed by various groups for a number of nuclear facilities that follow these guidelines.
Most relevant to the current project, the US Department of Energy recently submitted alicense application for the proposed Yucca Mountain geologic repository that uses
probabilistic methods and expert assessment to estimate probability and uncertainty
associated with volcanic disruption of that facility. Features of the methodologies common to
all of these analyses are:•

Alternative geological models of volcanism are considered and assessed. Any hazard
assessment should be consistent with geological models of the origin, ascent, andpotential eruption of magmas. Although there is broad agreement about the role ofplate tectonics in magmatism, details of these processes that are important on thescale of geological repositories are the topics of active research among the scientificcommunity.

•

Significant attention is devoted to verifying and assessing data used in volcanic
hazard assessments. This is especially important because widely available data, suchas geological maps, were not necessarily compiled for the purpose of assessment of
volcanic hazards. New methods are being developed to quantify the classification of
volcanoes, and hence their potential hazards, in a manner consistent with therequirements of volcanic hazard analysis.

•

Probabilistic methods are used to estimate the likelihood of volcanic events that may
disrupt the facility and to assess uncertainty in these estimates. The last decades
have seen tremendous improvement in the nature and capabilities of stochastic
models of a wide range of volcanic processes. These models have improved ourquantification of long-term volcanic hazards.

•

As in earthquake and fault hazard assessment, expert elicitation has become awidely used and accepted tool in volcanic hazard assessment. Formal expert
elicitation allows for a wide range of alternative geological models to be assessed byexperts with varying backgrounds and technical positions (often 10 or moreindividuals are involved in such assessments). Because volcanic hazard assessment
is a complex topic, such formal expert elicitations provide insights into theuncertainties extant in specific volcanic hazard assessments and sensitivities of theseassessments to specific assumptions or models.

A hazard assessment that includes these features can be readily used in performanceassessment to assess radiological risk, the primary metric used to assess site suitability.
Specifically, volcanic hazard assessment yields an estimate of the probability of volcanic
events. Once the probability of volcanic events is estimated, radiological risk can be
estimated using this value, coupled with analysis of the potential magnitudes of volcanic
eruptions, models of magma-repository interaction and models of transport of radionuclides
by volcanic processes. This risk-based approach is particularly important because of thetremendous variation in the scale of volcanic processes and their potential to result intransport of radionuclides into the biosphere (Figure 4.1).

57

ANALOGUES
Oshima-O-Shima volcano

GEOLOGICPROCESSES

MODELS

PA

Dike injection, conduit
development, radionuclidetransport in pyroclastic
flows and tephra

Megata Maars

Groundwater-magma
interaction, crater
development, radionuclidetransport in pyroclastic
flows and tephra

Hijiori Caldera

Development and eruptionof large volume magma
bodies, radionuclidetransport in pyroclastic
flows, tephra

Figure 4.1: New volcanoes formed in the Tohoku region during the last 50 Ka represent a wide range ofvolcanic events in terms of magnitude and potential effects of volcanism. These different magnitudes
and styles of volcanic activity should be reflected in performance assessment (PA) models of volcanism
for a proposed HLW repository site. Volcanic hazard assessment is primarily concerned withdetermination of the likelihood of such events.

This general framework, developed through numerous investigations by a large number ofgroups worldwide, is the basis for ITM methodology for assessment of the hazards and risks
associated with volcanism. ITM methodology development to date has focused on
assessment of the likelihood of volcanic events, using the Tohoku case study to elucidate therole of geological models, assessment of data uncertainty, and probabilistic methods in the
volcanic hazard assessment. Aspects of the Tohoku case study that are particularly relevant
to methodology development are discussed in the following.

4.1.1

Geological Models

Volcanoes in Tohoku are distributed along a volcanic arc and more diffusely distributed in a
back-arc region. In a classic view of subduction zone volcanism, the likelihood of new volcano
formation varies from near zero in the fore-arc, abruptly increases to a maximum across thevolcanic front, and decreases in the back-arc. An important feature of such a geologicalmodel is that along-arc variations in the probability of future volcanic events are minimal. The
Tohoku case study was used to evaluate this model and to consider more complex models of
volcanism in arc settings that have been developed during the last decade.

58

Figure 4.2: Arc magmatism (top) with zones of partial melting and their seismic tomography signature,
expressed in terms of S-wave velocity perturbations along the inclined low-velocity zone in the NE
Japan mantle wedge (bottom), correlated with the location of volcanoes (red triangles) after Hasegawa
and Nakajima (2004).

59Throughout the Tohoku case study, work by various Japanese groups (e.g. CRIEPI,JAMSTEC, Tohoku University, AIST) revealed correlations between volcanic phenomena andother types of geophysical and geological data. These data included seismic tomographic
anomalies, isostatic gravity anomalies, and topography (Figure 4.2). In the Tohoku region,
clustering of volcanoes occurs along the arc, and these clusters appear to extend into theback-arc, perhaps related to upwelling in the asthenospheric mantle in response to
subduction processes. These observations are the basis of the ‘hot finger’ model (Figure 4.3).
In this model, volcanoes are concentrated above finger-like zones of partial melting in theasthenosphere. Magma flux is high in these zones, under-plating the crust or perhaps leading
to sill development in the mid-crust. These density changes in the crust lead to isostatic
adjustments, manifest at the surface by up-warping of mountains in the back-arc of Tohoku.

Figure 4.3: Hot finger model according to Tamura et al. (2009). According to Kondo (2009), due to thestability of the plate tectonic regime and the continuity of the stress conditions for the past several
hundreds of thousands years, characteristics of future volcanism by analogy with Quaternary volcanism
are likely to remain similar at least for the next one hundred thousand years (up to c. 0.1 million years).

Consequently, from a conceptual point of view, the following pieces of evidence andhypotheses are integrated into ITM methodology: (1) the spatial distribution of Quaternary
volcanic events tends to be clustered rather than uniform along the arc, (2) the spatial
distribution of these events may be correlated with the geophysical signature of crustal andmantle structures, and (3) future volcanic events are more likely located in zones of past
activity and above hot regions within the mantle wedge. Our goal is not to prove or disprovespecific geological models, but to develop probability models that are sensitive to suchfeatures as clustering in volcano distributions and the presence of geophysical anomalies.
With such models, the affects of volcano clustering on estimates of the likelihood of volcanic
events may be assessed, and, if necessary, a weight assigned to specific models throughexpert elicitation.

4.1.2

Data Uncertainty

Studies of volcano distribution require a systematic and robust definition of a volcano indeveloping a database for analysis. Different definitions may result in different outcomes forstatistical modelling. If volcanoes in fact divide naturally into several groups, or alternatively
form a continuum, then these attributes need to be clarified.
The ITM volcanology group analysed volcanic hazard in the Tohoku region using data fromthe Catalog of Quaternary Volcanoes of Japan (Committee for Catalog of Quaternary
Volcanoes in Japan, 1999). This database is comprehensive, with information on location,age, volumes, rock types, dominant eruptive styles, morphology and miscellaneous otherfeatures, such as craters, intrusions, hot springs and fissures.

60

Unfortunately, definitions of volcanoes and volcanic activity vary across the region andnational databases are created with various objectives in mind, not only volcanic hazard
assessments. The Catalogue of Quaternary Volcanoes in Japan, which provides data on
volcanoes that have been active during the last c. 2 Ma (by the definition of the Quaternary inthe catalogue). Where any such database is applied for a specific purpose, it requires
interpretation. This database contains both temporal and spatial information, makes variableassumptions on how an individual volcano is defined and located, and may include biases
caused, for example, by artefacts of geological preservation potential.
Use of this catalogue uncovered some problems, mainly data inconsistencies. The internal
inconsistencies, which include different definitions of volcanoes and interpretations of volcanic
landscapes, led the group to form alternative datasets, more specifically, alternative volcano
definitions for volcanoes in the Tohoku region of NE Japan, which will be discussed in Section
4.2.

4.1.3

Probabilistic Assessment of Volcanic Hazard

Volcanic hazard assessment for a HLW repository needs to couple state-of-the-art geologicalmodels of volcanic processes with up-to-date geophysical data. Geophysical models will
always be imperfect. For example, the correlation between spatial density and timing of
volcanism are not well known. New details regarding the correctness of the ‘hot finger’ model
continue to emerge. Based on this evolving continuum of data and knowledge, it isinappropriate to develop a purely deterministic approach to volcanic hazard assessment, or torely too heavily on any one model or data set.
Statistics play an increasing part in the role of forecasting the probabilities of where and when
future hazards such as volcanic events may occur. Spatio-temporal models allow forecasts of
where volcanoes are most likely to form over designated timescales in specific regions.
Owing to the long performance periods required for geological repositories for radioactive
wastes, this type of forecasting is especially relevant in repository siting projects in regionsthat are clearly prone to volcanism.
As part of the Tohoku case study, the ITM volcanology group has explored three methods for
estimating the probability of new volcanism. One method investigates alternative volcano
datasets to generate 2D probability surfaces. This method uses a non-parametric Gaussian
kernel density function with varying isotropic bandwidths and tests the impact on probability.
Another method uses the same non-parametric kernel density function with anisotropic ‘bestfit’ bandwidths. These optimal bandwidths are estimated using three different optimization
algorithms that are based on the locations of the volcanic events used in the assessment. A
third method, the Cox Process model, was developed to include the use of multivariate
potential, which integrates geophysical data directly into the probabilistic analysis.
The latter two of these methods are compared at fourteen example locations, selected torepresent varying tectonic and volcanic environments within the case study region (describedin the Appendix to this report; Figure 4.4). This study analyses the variability in calculated
probabilities and attempts to quantify the variability and assess the uncertainty in the models
utilized. This approach is entirely new to volcanic hazard assessment and suggests thatdifferent locations, in different tectonic settings, are best assessed using different probability
models. As shown in the following, smooth probability models appear to work best in the
back-arc region, where volcanism is more diffuse. Models that are highly sensitive to clusters
appear to work best within the arc, where volcano clusters and gaps in volcanic activity areprevalent.

61

Figure 4.4: Distribution of the 14 example locations in the Tohoku region (numbered yellow boxes) withrespect to the location of Quaternary volcanoes (indicated by green circles); active faults are shown bydashed lines. Topography from Shuttle Radar Topography Mission (Jarvis et al., 2006).

The conclusion of the Tohoku case study for volcanism is a set of recommendations
developed for NUMO by the ITM Volcanism Group based on their findings. In the following,
the results are summarized and these recommendations are provided.

4.2

Alternative Data Sets and Cladistics

Defining a volcano or a volcanic event is not straightforward, but is a fundamental concern inthe creation or use of any volcano database. When using these data in hazard assessments,
it is particularly important to use volcano definitions that are robust and consistent. Within theCatalogue of Quaternary Volcanoes of Japan, volcanoes are commonly defined by theirmorphological descriptors, such as stratovolcano, maar, lava dome, caldera, or cinder cone.
Young volcanoes are easily mapped and catalogued; for example, the volcano location isdefined as the summit or main crater. However, difficulties arise when cataloguing older,eroded volcanic systems. In these cases, a more detailed study of near vent breccias or radialdikes is required to infer the appropriate location for the volcano (Connor and Hill, 1995).
Hone et al. (2007) discussed problems with dividing volcanoes in this database into
constituent “edifices” and minor inconsistencies in defining a volcano. The categorization of
volcanic features into volcanoes, and edifices belonging to a volcano, was based on the62interpretation of independent local experts; in particular, a volcano was defined as a single
coherent magmatic system (Prof. Sintaro Hayashi and Dr. Koji Umeda, personalcommunication) with edifices being part of that system. This process of human interpretationof volcanic landforms as objects for classification potentially introduces subjectivity andinconsistency, due to its qualitative nature and the involvement of many committee members(> 50) in the compilation of the database.
Classification schemes help to answer the question “what is a volcano?”, and are one methodof forming alternative datasets of volcanoes. Volcano classification schemes have previously
been largely based on either morphological characteristics or style of eruption. We havedeveloped a method employing the cladistics technique as an option of the classification
schemes. Cladistics is a systematic method of classification that groups entities on the basisof sharing similar characteristics in the most parsimonious manner. Cladistics is better known
for its use in breaking down all life into the ‘tree of life’, using a large dataset of characteristics
with no assumption being made about their relative significance. In this case cladistics takes
all the independent features of a volcano, compares each feature from each volcano withevery other volcano, then arranges them in the most parsimonious manner. Hone et al.
(2007) used all the characteristics of the 129 edifices and 59 volcanoes in the Tohoku dataset
to show that it can be divided into distinct clades. This type of analysis reduces human bias
in deciding on how volcanoes should be categorised and uses all the different parametersthat characterise each volcano. Volcano and edifice characteristics recorded in the databaseinclude attributes of volcano size, chemical composition, dominant eruptive products, volcano
morphology, dominant landforms, volcano age and eruptive history.

4.2.1

MethodData on the volcanic edifices were divided into 34 characters and 24 characters for volcanoes
based on the morphological features, size, style of activity, age and rock composition. Eachcharacter is coded which splits a characteristic into a small number of discrete bins, each ofwhich is given a code. The assembled data-matrix was run using software called PAUP*(Phylogenetic Analysis Using Parsimony, version 4.0; Swofford, 1998). The PAUP* code is
widely used by the biological community for systematic analysis.

4.2.2

Results

Table 4.1 provides a summary of the main features of the volcano groups formed fromcladistic analysis. In this exploratory analysis, different characteristics (e.g., the ages of
volcanoes) are included and different major groupings (1-4) identified through the cladistic
analysis. The result is a quantitative grouping of volcanoes with different characteristics. For
example, in analysis 4 the ages of volcanoes is not considered. Instead, three major groups
are identified based on composition, eruptive products, and morphological characteristics(groups 4a-4c). As these groups are not uniformly distributed along the arc (Figure 4.5) andare associated with different volcanic processes with potentially different consequences forrepository performance, their different distributions might impact hazard and risk assessments
at specific sites.

4.2.3

Clade-volcano cluster correlations

Spatial clusters of volcanic centres have been recognised in the Tohoku region by Tamura etal. (2002). The groups identified by cladistic analysis are distributed unevenly between theclusters, indicating a tendency for individual clusters to form similar kinds of volcanoes with
distinctive but coherent styles of volcanism. The cladistic analysis also shows that theclusters evolve with the character of the volcanism changing with time. Evolution ofcomposition from more mafic to more silicic volcanism is evident in clusters 4, 5 and 7. The
dominant magma composition affects eruptive styles, landforms and products. The cladistic
analysis appears to pick out these relationships because the time characteristics are linked tothe compositional controls. Some of the factors that might account for these differences may
include magma supply rate, longevity, crustal structure and magma rheology. The
classification results invite enquiry into the underlying dynamical causes of variability.

63

Table 4.1. Overall characteristics of the different cladistic groups (1A to 4C).
Age

Composition

Features

SlopesVolume3(km )1A

0-1.5 Mamixed

lava flows, pyroclastic rocks, lava
cones

steep

varied1B

Holocene

dacite to mixed

1C

0-0.5 Mabasalt to andesite

lava flows, pyroclastic rocks

steep to gentle

>100

1D

variedAndesitepyroclastic rocks, lava domes

gentlevaried2ANAbasalt to andesite

calderas, lava domes, lava flows,
pyroclastic rocks

gentle to steep

varied2BNAdacite to rhyolite

pyroclastic flow plateaus, lava domes

gentlevaried3Avariedmixed

lava domes, pyroclastic rocks

steep

10-100steep

10-100pyroclastic flow plateaus, lava domes,
steep to gentle
calderas

0.1-103B

0.5-1 Ma

dacite to andesite

calderas, lava domes, pyroclastic flow
plateaus

3C

>2 Ma

basaltic andesite to dacite

lava cones, lava domes, pyroclastic
flow plateaus

gentle

1-10

4A

NA

mixed

lava domes, lava flows, pyroclastic
rocks

steep

10-1004B

NAbasalt to andesite

lava cones

gentle1-10

andesite to rhyolite

pyroclastic flow plateau, lava flows,lava domes, calderas

4CNA

steep to gentle 10-100

Figure 4.5: Example of the application of cladistic analysis to the Tohoku volcano dataset. Three groups(4a-4c) are identified (see Table 4.1). These groups are distributed unevenly among the clusters
identified by Tamura et al. (2008). For example, group 4b volcanoes (predominantly smaller volume
basalt-andesite cones) are located in cluster 6. This suggests that styles, and hence consequences, of
volcanic activity are not uniformly distributed among the clusters. Such groupings, done in an
quantitative manner, can be considered in probabilistic assessments and risk estimates.

64

Our exploratory study indicates that cladistics has a promise as a method for classifying
volcanoes and potentially elucidating dynamic and evolutionary volcanic processes. A major
conclusion from this work is that the types of volcanoes, and hence style of volcanism, is notnecessarily uniform along the volcanic arc (e.g., Figure 4.5). Therefore, quantitative volcano
classification, for example using cladistics, should be used to estimate the probabilities offormation of specific types of volcanoes, with potentially different consequences for repository
performance. Traditional methods of volcano classification, such as those based onmorphology alone, are insufficiently precise for such an analysis. For specific sites, the affectof various classifications on estimated risk might then be assessed.

4.3

Probability Models

It is unlikely that any nuclear facility would be located in an area with a high probability of
future volcanism during its operational lifetime. However, geological repositories need toprovide isolation for tens of thousands of years and, over such timescales, the possibility of
future volcanism in areas that either have no previous record or are close to young volcanoes
needs to be evaluated. This is true for any volcanically active country or region in the world.
The constraints on siting a repository are that it should have an extremely low probability of
being directly impacted by future magma intrusion in the next 10 ka, with the probability ofsuch events over a 100,000 years timescale being an important aspect of the safetyassessment of potential sites.

4.3.1

Spatial Density

In the context of volcanic hazard assessments for repositories of radioactive wastes, thereason for estimating spatial density is to determine possible locations of future volcanic
events, or to estimate the probability of an event occurring at a specific location, given thatsuch events occur within the site region. Hazard assessments need to rely on robust regional
spatial density estimates of the likely occurrence of these future volcanic events. The
distribution of older volcanoes or volcanic vents provides some of the clearest information
about the probable location of future volcanoes and volcanic vents (e.g. Crowe et al., 1983;
Lutz and Gutmann, 1995; Connor and Hill, 1995; Jaquet and Carniel, 2006; Jaquet et al.,2008).
In the siting of nuclear facilities, however, estimating the spatial density of such hazardous
events has often proved contentious. Why the controversy? The underlying geologic
processes controlling the distribution of these events are complex and incompletely
understood. The frequency of such potentially catastrophic events is low, especially withinmany regions considered for nuclear facilities, so data used in these analyses are often
sparse. The selection of specific statistical models to estimate spatial density is often
subjective. These factors result in uncertainty.
Unfortunately, there is some ambiguity in the literature regarding the use of the terms density
and intensity. In the geosciences, variation in the number of events per unit area (say the
number of volcanic vents) is described using the term density. Intensity, in geoscience
contexts, often refers to the magnitudes of these events. The intensity of a volcanic eruption
can be characterized in terms of its total mass of eruptive products or related indices (Pyle,
2000).
Density and intensity are defined differently in spatial statistics. In this context, spatialintensity refers to the expected number of events per unit area defined at a point, s, a matrix
containing the x and y coordinates of the location of the point (Diggle, 1985; Diggle andMarron, 1988; Gatrell et al., 1996). Suppose there exists a set of events (e.g. volcano
locations) that occur within a given region, R . These events can be designated as
xn (n = 1,2,..., N ) ∈ R , where N is the total number of events, each consisting of the spatial
location, x and y of the event (possibly given in Easting and Northing coordinates, or, latitude
and longitude). One way we can create a model of spatial intensity from these events is toimagine they are realizations of a random variable, X , a function that describes the set of all65

possible realizations. For example, X might be the distribution of potential volcanoes, fromwhich a set of observed realizations (e.g. those found on a geologic map) are drawn. The
spatial intensity is formally written as (Gatrell et al., 1996)⎧ E ( X(ds)) ⎫⎬ds
ds → 0 ⎩⎭λ (s) = lim ⎨where E ( X(ds)) is the expected number of events that fall within a small area ds about the
point s (hence, if the location, s , is given as Easting and Northing with units of meters, thenthe units of λ (s) are m-2). At first glance it appears that the statistical definition of intensity isequivalent to the term density as commonly used in the geosciences. This is not quite true.
The geological processes that result in a given event distribution are incompletely known. We
can think of these geological processes as giving rise to a stochastic point process that
describes the relationship between the set of events and the geological processes that led totheir formation. As the stochastic point process is incompletely known, the value of thispotential local spatial intensity, λ (s) , is also unknown. That is, the observed distribution of
events is only one realization of the underlying process that gives rise to the set of potential
events. Our goals are to find an estimate of the spatial intensity, λˆ (s) , that approximates theunknown value of potential spatial intensity, λ (s) , and to understand the uncertainty in thisestimate.
In hazard assessments, there is a further requirement, that this information be used toforecast the spatial distribution of possible future events. Often we consider spatial intensity interms of the probable location of some future event, given that one occurs within our region ofinterest. This conditional probability can be estimated by:λˆ (s)d (s).
fˆ (s) =ˆ (s)d (s)λ∫RIntegrating fˆ (s) across the region of interest, R , gives unity, if R is sufficiently large. Since
all values of fˆ (s) within this region are greater than or equal to zero, this makes fˆ (s) aprobability density function and this function may be used in probabilistic hazard models.

fˆ (s) is referred to as one estimate of the spatial density, and one can consider the spatial
density per unit area in terms of conditional probability (e.g. given a volcanic event in theregion, what is the probability that the event will occur within some small area about the points ?). In addition, care is required in the selection of the region R , as external events locatedclose to the border may have a non-negligible contribution to spatial density. A practical
approach is to select R to be quite large compared to the region of specific interest (e.g. thesite).
How does one develop a best estimate of spatial density? In the real world, there is only one
realization of an underlying geologic process, the observed distribution of past events.
Unfortunately, geology is not conducive to repeating the experiment in a natural system. For agiven region there is just one geologic map of volcano distribution. Presumably, if there
existed a complete geophysical model for these events, we would use this information tobetter forecast the locations of future events. For example, if we knew the distribution of meltin the asthenosphere and lithosphere, and if we knew the state of the lithosphere throughwhich the magma rises, we might have a better sense of where volcanoes are most likely toform next. Currently, we lack such a complete geophysical perspective. Some data sets givean idea of where partial melting of the mantle might occur, for example seismic tomographic
models of "slowness" in the lithosphere and asthenosphere (e.g. Zhao, 2001). Other data,such as variations in gravity across a region (Connor et al., 2000; Parsons et al., 2006), show
some correlation with the existing distribution of volcanoes in some circumstances, but the66

mechanisms relating gravity anomalies to the origin of magmas are not completely
understood. As a result, these types of data have been used to support estimates of spatial
density (e.g. Connor et al., 2000; Martin et al., 2004; and the Cox Process model described in
the following), but no model has yet been proposed that does not rely principally on thespatial distribution of past events.
The reliance on the distribution of past events implies that these realizations are
representations of some underlying random variable, X , that will govern the distribution ofpotential events in the future. This assumption immediately raises a fundamental question.
Which are the past events that should be used to develop the spatial intensity estimate, λˆ (s) ,and density, fˆ (s) ? Event datasets used to estimate the spatial density of future events need
to be consistent with several features of geological processes.
First, any spatial intensity function for a geologic process must change with time. On timescales of tens of millions of years, plate boundaries change, volcanic arcs wax, wane, andmigrate, and major fault systems reorganize. In very long term probabilistic hazard
assessments for high-level waste repositories, which may have 106 years of performance
periods, these factors have to be considered in weighing the validity of using specific data indeveloping spatial intensity models. For processes like volcanism, where a geologic record of
past events usually persists for tens of millions of years, consideration needs to be given towhich events best represent the distribution of future volcanism (e.g., Kondo, 2008). For
example, the distribution of Miocene volcanoes in a given area might be much less relevant
than the distribution of Pliocene and Quaternary volcanoes. Thus, in order to develop an
estimate of the spatial intensity, a model of the geologic evolution of the system is required.
This geological model is used to justify the inclusion of some geological features in the event
dataset, and the exclusion of others.
Second, it is necessary to assess the completeness of the geologic record. Volcanic eventsmight be missed in initial geological investigations, as volcanic vents might be buried in
sediment or otherwise obscured (e.g. Connor et al., 1997; Wetmore et al., 2008).
Third, geological events, even when they are all identified, may be so rare as to present an
incomplete picture of the underlying process. Consider the formation of a new volcano as asingle event, x n , one realization of the random variable, X . If, for example, X can becharacterized by a uniform random distribution, then it is likely that the observed set of
realizations will have a spatially random distribution within the region of interest, R . However,the underlying density usually has additional structure, causing independent realizations tocluster. For example, volcanoes cluster above zones of partial melting in the mantle. For
random variables with a great deal of statistical structure, such as many modes in spatialintensity, a great number of events might be required to identify the statistical structure of therandom variable.
Fourth, it is critical to ascertain which geologic features are actually independent events. The
true statistical structure of the random variable, X , might be obscured if some events
included in the event dataset are not independent. Volcanoes are complex geologic structures.
The spatial distribution of polygenetic volcanoes reflects processes of magma generation andrise through the crust. The distribution of small vents (sometimes referred to as parasitic oradventive cones) does not necessarily reflect the distribution of polygenetic volcanoes, so aspatial intensity estimate that includes all vents as events would not correctly model theunderlying random variable. Furthermore, in monogenetic volcanic fields alignments of
volcanic cones develop in response to single magmatic events, episodes of magma risethrough the shallow crust. This is because single igneous dikes ascending through the crustmight form segments and rotate within the shallow crust, each segment feeding a separate
vent and each building a volcanic cone. If the goal of analysis is to forecast the distribution of
future magmatic events, each of which might produce more than one monogenetic volcano,geological data must be gathered and volcanoes formed by the same magmatic event must
be somehow grouped as single events.

67Independence of events is not necessarily easy to determine. Rather than simply counting
volcanoes on a geologic map, one must make a geologic assessment of the independence of
these data. For volcanoes, this is generally accomplished through detailed analyses ofradiometric age determinations, stratigraphic correlations, and related geologic data. Often,
even detailed analyses do not resolve whether or not specific features should be grouped as
single events or treated as separate, independent events. Consequently, a major task in
preparing a spatial intensity estimate for a nuclear facility is defining the dataset of events to
be used. Certainly a major expense in site characterization is data gathering to support
interpretation of geological features as events.
Hazard estimation requires the estimation of spatial density models. As described previously,
hazard assessments often consider alternative event datasets and account for the affect ofthese varying datasets on spatial density estimates. The ITM volcanology group uses two
related probability models, the non-homogeneous Poisson process with specialized kernel
methods and the Cox process with multivariate potential, for characterizing the potential of
volcanism. The Cox process assumes the potential to be a random intensity function in orderto account for its uncertainty with respect to hazard modelling. Both statistical models aredesigned to reflect basic features of volcanism in Japan, specifically the clustering of volcanic
centres. As the models are different in their underlying assumptions, the estimated hazard
from application of these models will differ. This difference gives an indication of epistemic
uncertainty in hazard assessments. In the following we summarize the specifics of these twomodels in order to characterize the potential of volcanism.

4.3.2

Estimating Spatial Density with Kernel Methods

Within the framework of the non-homogeneous Poisson process, the probability of new
volcanic events occurring within a small domain of the region of interest is given by:

P { N i = n} = e− λˆs , i A i(λˆ A )
s, i

n

i

n!

where Ni is the number of volcanic events within the domain Ai partitioning the region ofinterest. The characterization of the deterministic potential is achieved by kernel estimation
techniques, widely used to model spatial distributions of volcanoes. In this technique, spatial
variation in

λ̂s is cast as a function of distance to nearest-neighbour volcanoes and asmoothing parameter, h, known as the kernel bandwidth. The parameter h may be a matrix H,
which allows the smoothing to be asymmetric and sensitive to geologic features. These
analyses may be weighted by geophysical data to account for incomplete correlation between
geophysical data and spatial density of volcanoes. Uncertainty in the models is assessment
using a smoothed bootstrap technique and Monte Carlo simulation (details are given in
Connor and Connor, 2008)

4.3.3

Cox process with multivariate potential

The concept of the potential of volcanism (Jaquet et al., 2008) is defined as the propensity ofa given region to be affected in the future by volcanic events. Using this concept allows for thedescription of event distribution (e.g. clustering effects) as well as for the integration of
additional relevant data (e.g., geophysical data).
Within the framework of the Cox process with multivariate potential, the probability of new
volcanic events occurring within a small domain of the region of interest is given by:⎧
Zn ⎫P { N i = n} = E ⎨e − Z i i ⎬
n! ⎭
⎩

68

The volcanic region of interest is partitioned into small domains, where the random potentialof domain Ai is denoted by Zi. In comparison to the non-homogeneous Poisson process: (1)the potential corresponding to the mean number of volcanic events in domain Ai is interpreted
as a realization of a random intensity function and (2) the numbers of volcanic events within
disjoint domains are no longer independent, due to the structured behaviour of the random
potential modelling the observed geological patterns. Among the contributing factors to thepotential, the location of zones of partial melting at depth is of interest for hazard assessment.
Such information can be provided by geophysics, which in the model is also interpreted as arealization of a random function. The dependence relationship between the potential and thegeophysical information is modelled under the assumption that their respective transforms in
Gaussian space are related by a regression formula involving a statistical coefficient ofcorrelation (details are given in Jaquet and Lantuéjoul, 2008). In the developed model, thepotential of volcanism becomes multivariable, since it presents dependencies with past
volcanic activity as well as with geophysical data.
The estimation of volcanic hazard is performed by simulating the distribution of volcanic
events likely to occur during a certain period of time in the future within the region of interest.
In addition, the simulation has to deliver volcanic events that are more likely to be located inzones of past activity. Therefore, the simulation requires to be conditioned to the number ofpast volcanic events known in each domain ai . The idea is to simulate the potential of
volcanism conditioned on the number of volcanic events and on the geophysical information
known for the analyzed region. This simulated potential is representative only of the period oftime from which all data originate. For the simulation of volcanic events, occurring in a future
period of time, the future potential is required. Provided that the potential varies very slowly
through time, there is no inconvenience to assume that past and future potentials areproportional. The details of the proposed algorithm for the conditional simulation of thepotential are given in Jaquet and Lantuéjoul (2008).

Figure 4.6. Flow chart of the Cox Process methodology for the estimation of volcanic hazard.

69

The conditional simulation algorithm allows the estimation of volcanic hazard for each domainof the region of interest during the period of time considered. A Monte Carlo approach (cf.
Figure 4.6) is performed using several thousand simulations in order to derive stablefprobability estimates for the future volcanic events, N i :

P

{N

f
i

≥1

}≈ K1 ∑ 1 ( n
K sim

k

sim

f
i

≥1

)

k =1

(

)

where K sim is the total number of simulations and 1k ni ≥ 1 equals 1 when the kfthsimulation assigns the domain ai one or more volcanic events, and 0 otherwise.

4.4

Application of Probability Models in Tohoku

We apply these models in the Tohoku region in order to assess their performance withrespect to:•Modelling the distribution of volcanoes (clusters)•

Consistency with geological models of volcano distribution in the Tohoku region (e.g.,
the classic volcanic arc model and the hot finger model of volcanism, which results involcano clustering)•

Success in forecasting the potential of volcanism evaluated through the blind test to
identify where volcanoes formed during the last 50,000 years based only on the
distribution of older Quaternary volcanoes.

The goal of this models development is not to choose among various models but tounderstand the range of hazard estimates that emerge from application of various wellconstructed statistical models.
The issue of concern is the estimation of the probability for the formation of new polygenetic
volcanoes over a possible performance period of 0.1 million years. All Quaternary polygenetic
volcanoes (those active within the last 1.8 million years) of the Tohoku region are considered
as representing potentially active volcanic areas for the future. A volcanic edifice corresponds
to the volcanic event considered for the Tohoku case study; i.e. a polygenetic volcano iscomposed of one or more edifices. Polygenetic volcanoes are characterized by ageomorphology and geologic structures that are created by many episodes of eruptive activity
likely to affect broad areas. The formation of new polygenetic volcanoes can take place at
locations up to tens of kilometres away from sites of previous eruptions. Some of the
volcanoes of the Tohoku region remain active today. We illustrate our methodology using a
subset of the Database of Quaternary Volcanoes of Japan. Subsets of these data, such as
those identified through cladistic analysis, are not included in the following only because thecurrent research has focussed on model development and comparison of models with astandard data set. In an actual hazard analysis for a specific site, numerous alternative data
sets should be considered.
The potential of volcanism for the Tohoku region is believed to assimilate geological data
related to Quaternary volcanic events as well as seismic tomographic models indicatingzones of partial melting at depth. The latter consists of P- and S-wave velocity perturbations
along the inclined low-velocity zone in the mantle wedge of NE Japan (Hasegawa andNakajima, 2004). For a specific site, other geological and geophysical data sets should beconsidered, depending on geological models. The same procedures for correlating thesealternative geological and geophysical datasets may be used.

70Probabilities are estimated using: (a) spatial densities from the various statistical models, (b)homogeneous recurrence rates for Quaternary volcanism, (c) areas that reflect site vicinity,
within which the effects of potential volcanism are considered to be deleterious to repository
performance (say 25 km2), and (d) performance periods (say 100,000 yr). Volcano event
magnitude is not considered (e.g., large events, such as caldera formation, may affect region> 25 km2 but this is not considered in the example model). The average Quaternary
recurrence rate of volcano formation in the map region is 103 events / 1.7Ma = 6 x 10-5 events
per year. Thus the probability of an event is calculated using a non-homogeneous Poissonmodel for low probability events:

P { N i ≥ 1 |Δt , A} = 1 − e− λˆ s ,i λˆ t Δ t A

where Ni is the number of events expected to occur in time period, say, Δt = 100,000 years
within a domain A= 25 km2, with estimated spatial density of volcanic events, λ̂s andestimated temporal recurrence rate of events, λ̂t = 6 x 10-5 events per year.

4.4.1

Kernel Density Estimation

Non-parametric kernel density functions using anisotropic bandwidths estimated by
optimization algorithms were used to assess spatial density of volcanism in the Tohoku region.
These methods include the Least Square Covariance (LSCV), Smoothed Covariance (SCV),
and Smoothed Asymptotic Maximum Squared Error (SAMSE) methods. Althoughmathematically complex, these methods are all designed to achieve the basic goal of using
the volcano distribution data in the most efficient way possible to estimate spatial density, andhence probability.
The methods all use the distance between volcanoes to estimate spatial density. The SAMSE
and SCV methods use a smoothed representation of these distances, whereas the LSCV
method uses the unsmoothed distances. This means that the SAMSE and SCV methods tend
to produce a smoother map of spatial density than the LSCV method, which is highly sensitive
to clustering.
Figure 4.7 shows the results of applying these algorithms to estimate spatial density of
volcanism in the Tohoku region. The SAMSE and SCV methods produce smooth models ofprobability, emphasizing the lack of volcanism in the fore-arc, the high concentration of
volcanism at the volcanic front, and the lower concentration of volcanism in the back-arc
region. These models are most consistent with ‘classic’ subduction zone and back-arc basinmodels, which emphasize across arc variations. Note that these smooth models are notparticularly sensitive to clustering of volcanism within the arc, or along arc variations in spatial
density.
The LSCV probability models (Figure 4.7), like the Cox Process models discussed in thefollowing section, emphasize volcano clustering within the arc and behind the arc. Thesemodels are most consistent with the ‘hot finger’ model, which predicts gaps in volcanism
along the arc. These models are also most consistent with the distribution of mantle seismic
tomographic anomalies. Different kernel density estimates of spatial density are sensitive todifferent aspects of the tectonic setting of volcanism in Tohoku.

71

Figure 4.7: Contour plots show the probability of volcanism based on spatial density estimates madeusing a Gaussian kernel function and three mathematically derived ‘optimal’ bandwidths, identified as
LSCV, SAMSE and SCV. Darker colours indicate areas of higher probability. The boxed numbers
indicate the locations of the 14 example locations (see Figure 4.4). Note that the SAMSE and SCV
models yield smooth results, consistent with classic models of subduction zones. The LSCV model, likethe Cox Process model described in the following, is more highly clustered, more consistent with the ‘hot
finger’ model of volcanism. On all maps, darker colours indicate higher probabilities of volcanic events.
Probabilities estimated for a 25 km2 area and 100,000 yr period.

4.4.2

Cox Process Model

The simulation of the potential using a Cox Process model requires the estimation of thefollowing model parameters: (1) the Gaussian transforms for the potential and the seismic
data, (2) the correlation coefficient and (3) the variograms describing the randomly structured
spatial patterns for the potential and the seismic. The details of the estimation procedure foreach of these parameters are given in Jaquet and Lantuéjoul (2008).
With the Cox model, simulations of the number of future volcanic events were carried out overa period of 0.1 million years for the Tohoku region (Figure 4.8). A Monte Carlo approach was
performed using 10 000 simulations in order to obtain stable probability estimates. The
likelihood of future volcanic events was displayed in form of a hazard map for the period ofinterest related to the siting for a HLW repository (Figure 4.9).

72

Figure 4.8: Three Cox simulations with a multivariate potential of volcanism. The simulated edifices(large white triangles) are likely to be located in zones with past activity as well as in zones with seismic
anomalies (i.e. low dVs values derived from seismic tomography models). Yellow to red colours indicate
lower seismic wave velocities, indicating increased potential for the presence of partial melt in the
mantle.

Figure 4.9: Probability map of new volcano formation in the Tohoku region and for the next 100,000 yr(domain scale is 10 km by 10 km). Darker orange colours indicate a higher spatial density, and therefore
higher probability of a volcanic event.

73

4.5

Assessment of Probability of Volcanic Disruption at 14 Example
Locations

The probability models developed by the ITM Volcanism Group were tested using 14 example
locations selected by ITM in the Tohoku region (Figure 4.4). These 14 locations were selected
merely to assess probability model performance; they have no programmatic significance for
NUMO. Locations were selected within fore-arc, arc, and back-arc regions, within volcano
clusters, 15 km from volcano centres, and farther from active Quaternary volcanoes. Someare located with “gap” regions within the volcanic arc. Additional details about individual
locations are available in the Appendix to this report.
Maps of the variation in spatial density across the region of the 14 locations are shown in
Figures 4.10 and 4.11 for SAMSE and LSCV estimates of spatial density, using the
nonparametric kernel density estimate. Table 4.2 compares the probabilities estimated at the
14 locations using SAMSE and LSCV kernel models, and the Cox Process model with amultivariate potential.

Figure 4.10: A smooth (SAMSE) model of volcano spatial density calculated for the Tohoku region andshowing the 14 example locations. This model best reflects a classic subduction zone model, with little
variation in spatial density along the volcanic arc. Spatial density is shown as the probability of a new
volcano forming within a 25 km2 area during 100,000 years. Darker orange colours indicate a higher
spatial density, and therefore higher probability of a volcanic event.

74

Figure 4.11: A clustered (LSCV) model of volcano spatial density calculated for the Tohoku region andshowing the 14 example locations. This model best reflects a ‘hot finger’ model, with substantial
variation in spatial density along the volcanic arc. Spatial density is shown as the probability of a new
2
volcano forming within a 25 km area during 100,000 yr. Darker orange colours indicate a higher spatial
density, and therefore higher probability of a volcanic event.

At the 14 locations, the three different models of spatial density result in probabilities that varyby approximately one order of magnitude. The Cox Process and LSCV models yield quite
similar results at most locations. An exception is location 14, located relatively far from active
Quaternary volcanoes but located within a zone of relatively low seismic velocities, possibly
indicative of the presence of partial melt in the mantle. For a potential site at such a location, it
is clear that further investigation of the origin and nature of these geophysical anomalies
would be necessary. As the SAMSE model yields a more smooth result, probabilities resulting
from the SAMSE spatial density estimate tend to be lower in the most active volcanic clusters,
and higher in gap regions within the arc between clusters, and in the back-arc region. The
smoother nature of the SAMSE solution is most evident at location 9, located in the fore-arc.
Because the SAMSE solution is smooth, the probability is arguably too high at this site. On
the other hand, the Cox Process model and the LSCV model yield low probabilities
throughout much of the back-arc region. These spatial density models tend to emphasize
clusters of individual volcanoes or pairs of volcanoes (Figure 4.11). These models may place
too much emphasis on clustering in the back-arc, whereas the SAMSE spatial density model
produces a smooth solution in the back-arc that appears to be more consistent with observed
volcano distribution (Figure 4.10).

75

Table 4.2: Comparison of the spatial density estimates for the 14 locations (see Figure 4.4) and using
different methods to estimate spatial density. In all cases, an area of 25 km2, a period of 100,000 yr, anda temporally homogeneous recurrence rate of 6.5 x 10-5 events per yr is used. Thus, the only factor
varying among the analyses is the estimate of spatial density of volcanism.
Example Location

Probability and Method

Remark

Cox

LSCV

SAMSE

-5

1.7·10-3

1

Within gap region

3.7·10-4

2

Within gap region next to fault

6.5·10-4

7.9·10-4

2.7·10-33Within Sengan cluster

2.6·10-2

1.1·10-2

3.1·10-3

4

15 km to west of cluster

7.7·10-4

9.7·10-4

9.9·10-4

515 km from center of cluster

8.1·10-3

5.9·10-3

2.2·10-3-3-3

1.4· 10-3-101.6·10-4

1.6·10-4

7.2·10-4

6.1·106

15 km to east of cluster

3.5·103.1·107

Sea of Japan site: outside modelled domain for Cox

NA

3.5·108Noshiro site

3.7·109Hachinohe site: outside modelled domain for Cox

NA-4

2.9·10-88.3·10-5

1.1·10-3

1.2·10-3

7.3·10-3

2.3·10-3

01.6·1010

Odate site

3.5·10-4

11

Iwate site

3.6·10-2

12Miyako site: outside model domain for Cox

NA13

gap site

3.5·10-4

1.6·10-3

2.0·10-3

14

outside 15 km exclusion zone

1.2·10-4

6.8·10-6

9.5·10-4-11This result indicates that it is useful to use a range of probability models to evaluate volcanic
hazard at a specific site, and to assess the performance of these models through comparison.
The relative weight given to specific models might vary depending on the nature of the site
region. In this case, analysis suggests that spatial density models that are sensitive toclustering are most useful in arc and fore-arc locations, and smoother models might be more
appropriate in back-arc locations.

4.6

Sensitivity to Probability Models

The ITM methodology involves a risk-based analysis of volcanism, with the goal of
integrating volcanic hazards into Performance Assessment (PA). Probability models shouldbe evaluated in terms of statistical confidence, and in terms of consistency with a range ofvolcanological concepts and models. This is apparent at the transition from development ofgeological models to their simplification and implementation in performance assessment (PA).
Integration into PA requires four steps:(1) Assessment of probabilistic forecasts of volcanic events(2) Assessment of the potential magnitudes of future volcanic events(3) Modelling of magma interaction with repository engineering and geologic systems
and potential dispersion of radionuclides into the biosphere by volcanic processes(4) Abstraction of these models in PA to estimate risk
In the Tohoku case study we have begun to consider points (1) and (2). One way to assess
probability models is to consider a subset of the data, and investigate the performance of themodel using the remaining data to forecast the occurrence of this subset. We select a subset
from the entire Quaternary database consisting of volcanoes thought to have formed within
the last 50 Ka. These volcanoes are removed from the Quaternary dataset and hazard
models re-calculated. Then the performance of the models is evaluated by comparing
forecast probabilities based on these models with the positions of volcanoes formed during
the last 50 ka. One would hope that probability models would do a much better job thanrandom chance of forecasting the positions of volcanoes in the subset.

76

In the greater region (Sapporo to central Honshu), 16 volcanoes or volcano groups formed in
the last 50,000 years, based on the oldest radiometric dates available for these volcanoes(AIST Quaternary database) (Table 4.3).
Table 4.3: Sixteen volcanoes formed in the greater Tohoku region during the last 50 Ka. Finger number
refers to the volcano cluster, as grouped by Tamura et al. (2002). Distance to nearest-neighbour older (>50 Ka) Quaternary volcano is given in km. Most volcanoes formed within 15 km of their nearest
neighbour (lighter-colour), but four formed >15 km from their nearest neighbour (darker-colour).
Finger Number

Volcano Name

Near-neighbour Distance (km)10

Toya-nakajima

0.0

6

Megata2.3

1

Nantai-san

3.6

10

Usu-zan

6.4

9

E-san

8.0

10

Eniwa-dake

8.4

4

Hijiori Caldera

9.8

5

Naruko Volcano Group

12.1

9

Zenikame Caldera

12.5

10

Shikotsu Caldera

13.4

6

Kampu-zan

13.9

9

Koma-ga-take

13.9

10

Fuppushi-dake

18.6

10

Tarumai-san

20.2

9

Nigorigawa Caldera

23.1

8

Oshima-Oshima

40.5

This rate of activity suggests a late Quaternary (< 50 ka) rate of volcano formation of
approximately 15/50000 = 3 x 10-4 per year (Toya-Nakajima volcano is excluded from theanalysis because it formed within an older caldera and has 0 km distance to nearestneighbour, although it is designated as a separate volcano in the AIST Quaternary volcano
database). This is apparently a much higher rate of activity of new volcano formation than in
the Quaternary as a whole. It is possible that that this rate of new volcano formation simply
reflects radiometric age determinations. Radiocarbon age determinations best estimate ages
to approximately 50 Ka, whereas K-Ar age determinations tend to bias estimates to older
ages. Therefore the higher 50 Ka recurrence rate may simply reflect sampling bias resulting
from radiometric age determination techniques.
It is also likely that some of the new volcanoes reported in Table 4.3 are actually part of a preexisting volcano plumbing system. Presumably, the closer the “new” volcano is to an older
Quaternary volcano, the more likely it is to be associated with an existing magma plumbing
system. If one uses the 15 km distance from older Quaternary volcanoes as a simple criterion,
then 4 new volcanoes formed in the region during the last 50 Ka. Applying this criterion yields
a rate of 4/50000 = 8 x 10-5 events per year, a value much closer to the average Quaternary
recurrence rate of 6.5 x 10-5 events per year. For the purposes of volcanic hazard assessment
for a specific proposed site, the lower average recurrence rate makes more sense to apply,
given that NUMO will not consider potential sites closer than 15 km to an existing Quaternary
volcano. This suggests that the higher rate may reflect some eruptions from pre-existing
volcanic systems, and not the formation of “new” volcanoes.
Figures 4.12 and 4.13 illustrate the success of two of the spatial density models (SAMSE and
LSCV) in forecasting the future positions of volcanoes. Most volcanoes form in areas ofrelatively high spatial density. It is most interesting to note that the LSCV (and Cox Process)77models perform better than the smoother probability models (SAMSE and SCV) in forecasting
probabilities where new volcanoes form within and close to the volcanic arc. That is, because
volcanoes cluster within the arc and new volcanoes tend to form within these clusters, higher
probabilities are forecast by the LSCV model where the volcanoes actually formed.
Conversely, the smooth probability models tend to overestimate probabilities of new volcano
formation in gap regions, between volcano clusters within the arc. In contrast, the smootherprobability models appear to perform better in the back-arc region, where volcanism is lessdensely distributed and the tendency toward clustering is less prevalent. Although the “hotfinger” model suggests that clustering occurs in the back-arc, it is less obvious in the
distribution of volcanoes than might be hoped. Consequently the smoother probability models
are less likely to underestimate future volcano distribution in the back-arc, at least based on
the formation of volcanoes during the last 50 Ka.
It is equally important to consider the magnitude of volcanic activity associated with new
volcano formation. Volcanoes formed in the last 50 Ka include a very large range of eruptive
magnitudes and potential effects. These volcanoes include: composite volcanoes (e.g., 1790
AD Oshima-O-Shima composite volcano – closely associated with older pre-Quaternary(?)
volcanoes Higashi-yama and Nishi-yama, which are not listed in the AIST Quaternary volcano
database); maars and volcano groups (e.g., The ~30 Ka Megata Volcano Group comprises
maars including Ichinomegata, Ninomegata and Sannomegata), and calderas (e.g., The ~12
Ka Hijiori Caldera, 2.5 km diameter). Thus, there is tremendous variation in the magnitude of
volcanic events during the last 50 Ka. Hazard models should encompass this range oferuption magnitudes and their subsurface effects. Event magnitude may be included inhazard assessments through study of analogue volcanoes, numerical simulation of volcanic
processes (e.g. tephra fallout and transport of radionuclides), and abstraction in PA with
Monte Carlo methods (e.g., Figure 4.1).
In summary, major findings of the ITM Volcanism Group with respect to tests of probabilistic
models and the magnitudes of potential events, based on the Tohoku case study include:•

All of the probability models tested are better at forecasting the distribution of
events over the last 50 ka than completely spatially random models (Figures 4.12
and 4.13).

•

The cluster models (e.g. LSCV, Cox process) perform better than the smoothmodels (SCV, SAMSE) overall for forecasting events that occurred during the last
50 Ka. This suggests these models might be given higher weight in forecasting
future events. An exception is forecasting volcanism in back-arc areas, where thesmooth models appear to perform better.

•

The formation of one volcano, Oshima-O-Shima, was not well forecast by thesemodels, although there are older, pre-Quaternary (?) volcanoes in this area. In
site–specific investigations, pre-Quaternary volcanoes should be considered in
the volcanic hazard analysis.

•

Very large magnitude volcanic events occurred in the region during the last 50 Ka(e.g. formation of 4 calderas). Hazard thresholds should take such large
magnitude events into account, and a full range of potential events should beconsidered in volcanic hazard assessment and PA. It may be necessary at
specific sites to quantify volcano classification, for instance using a cladisticclassification, and analyze the probability of volcanoes associated with specificeruption styles.

78

Figure 4.12: Success of the probability model in forecasting future volcanic events, using a smooth(SAMSE) model of volcano spatial density calculated for the Tohoku region. Older Quaternary
volcanoes (green circles) are used to estimate spatial density. New volcanoes (thought to have formed <50 ka) are indicated by red circles. This type of smooth model performs better for locations in the backarc. Spatial density is shown as the probability of a new volcano forming within a 25 km2 area during
100,000 years. Darker orange colours indicate a higher probability of new volcano formation estimatedfrom the statistical model.

79Figure 4.13: Success of the probability model in forecasting future volcanic events, using a clustered(LSCV) model of volcano spatial density calculated for the Tohoku region. Older Quaternary volcanoes(green circles) are used to estimate spatial density. New volcanoes (thought to have formed < 50 ka) areindicated by red circles. This type of clustered model performs better in forecasting future volcano
events within and near the volcanic arc, a region characterized by clustering and gaps in volcano
distribution. Spatial density is shown as the probability of a new volcano forming within a 25 km2 area
during 100,000 years. Darker orange colours indicate a higher probability of new volcano formation
estimated from the statistical model.

80

5 Conclusions and Recommendations
The Tohoku Case Study is the first test of the ITM methodology and, indeed, was used as atest bed to develop the methodology. Consequently the results should be seen as a first
demonstration of the types of information that a probabilistic approach to tectonism canproduce and, we hope, a useful illustration of how this information might be used by NUMO insiting decision-making and safety evaluation.
Large parts of the methodology are novel in the field of tectonic hazard assessment and thetotal approach has been established and tested over a relatively short period of about three
years. The ITM project is currently extending the methodology to another Case Study area,
after which there will be a broader experience base upon which to judge the usefulness of thetechniques and approaches developed.
Overall, the ITM team believes the Tohoku Case Study to have proved extremely successful.
It has shown that quantitative methods can be used not only to discriminate between
locations scattered within a given region in Japan, but also to identify 'quiet' locations where
susceptibility to tectonic impacts in the next 10,000 or even 100,000 years will be quite low.
Whilst this may have been empirically expected at the beginning of the study (e.g. from thepresence of 'gaps' in volcano distribution), the ITM methodology provides these observations
with objective and quantitative support.
In the following sections we highlight some specific observations and recommendations
arising from the rock deformation and volcanics teams.

5.1

Rock Deformation

The rock deformation probabilistic methodology, using various tectonic datasets of the
Tohoku region, has generated significant new knowledge and demonstrated an approach thathas not been achieved previously. While the method is new, it is built from adaptation of wellproven and widely used seismic hazard methodology.
The general agreement between the strain rates from GPS, which is a measure ofcontemporary strain accumulation, and those from seismicity and surface deformation,
derived from a long-term geological record of strain release, is particularly interesting because
it may suggest that the GPS strain is largely elastic and will be released in future deformation
events. We need to assess whether this result holds true in other tectonic domains,particularly where crustal extension occurs in association with volcanism.
In the Tohoku case study we have undertaken some preliminary uncertainty analysis. At thisstage of the case study we are aware that there are numerous further datasets and alternate
geological models that would need to be considered if volunteers were to come forward in
Tohoku. We have been assisted by Japanese experts in a trial expert elicitation process ofalternative models and parameter weighting. However, the opportunity has been limited, andseveral aspects of data and analysis presented here have not been through even thispreliminary expert elicitation process.
Although we have assumed rather arbitrary values for unfavourable rock deformation
characteristics in relation to example location characterisation, the procedures we havedeveloped do appear to offer an excellent basis for comparisons and for analysing particular
aspects of rock deformation hazard. The procedures and methodologies developed with the
Tohoku case study must be considered to be ‘alive’ in the sense that as new or improved
information becomes available, then iterations of the procedure must be undertaken to refine
the hazard assessments and associated uncertainties.

5.2

Volcanism

Throughout the Tohoku case study there were questions and concerns about databaseconsistency and the meaning of particular database elements. In the earliest stages of81assessing volunteer sites, it will be necessary to develop a volcano database specificallydesigned for volcanic hazard assessment. This database can draw on national databases(developed for a variety of reasons), volcanological literature, and maps. We recommend that
this database be developed by volcanologists familiar with the site region(s). Furthermore, it isrecommended that hazard analyses be made with alternative databases, and that the results
obtained using alternative databases be compared and assessed.
In addition to geological information, a variety of geophysical information should beconsidered in volcanic hazard model development for specific sites. This information includes
geophysical data sets that are available in the literature, such as seismic tomography, gravity,and perhaps electromagnetic soundings. This information includes geophysical models of arc
volcanism, back-arc volcanism, and volcano clustering.
Probability models should be used in the volcanic hazard assessment. Probability models are
standard internationally and have been effectively applied in the Tohoku case study. A variety
of probability methods have been developed as part of the ITM case study in Tohoku. It is
recommended that this range of methods be used in site evaluation in order to evaluate
uncertainty associated with the application of specific statistical models. It is noted that in the
Tohoku case study, the utility of specific probability models varied depending on tectonic
setting (arc or back-arc site regions).
The regional recurrence rate of volcanism appears to be a major source of uncertainty,
primarily because of difficulty in determining the age of formation of specific volcanoes. It is
recommended that conservative estimates of the regional recurrence rate of volcanism beused in site evaluation. It is recommended that the PIA stages of site assessment include a
program of radiometric age determinations to supplement geological mapping in the site
region, in order to acquire better estimates of recurrence rate.
In the Tohoku region, most volcanoes that formed within the last 50,000 years formed lessthan 15 km of an older Quaternary volcano. This provides direct support for NUMO’s 15 km
exclusion zone as a method of minimizing hazard and risk – not only from eruptions ofexisting volcanoes, but for the formation of new volcanic centres as well. Nevertheless, somevolcanoes formed during the last 50,000 years more than 15 km from older volcanic centres.
Based on the Tohoku study, it appears that detailed volcanic hazard assessment will benecessary for any site located in the arc or back-arc region.
It appears that at least in one location in the Tohoku region, a volcano formed during the last
50,000 years that was outside exclusion zones, but very close to pre-Quaternary (?)
volcanoes. It is recommended that such cases be considered during site evaluation, and thatthe hazard assessments and models include Pliocene volcanoes.
Data used in this case study are basically restricted to the Catalogue of Quaternary
Volcanoes in Japan because it was the most reliable and comprehensive database available
at the beginning of the study. It is obvious from the above discussion that incorporating data
older than the Quaternary to the dataset is essential for the assessment of future volcanism.
Analysis of volcanism in the Tohoku region, particularly over the last 50,000 years, indicatesthat a wide range in the magnitude and nature of volcanic events must be considered inhazard analyses. Hazard models should rely on a range of eruption scenarios and numericalmodels of volcanism to estimate risk. These models may evolve and become moresophisticated in PA assessments, if the hazard warrants such attention.
Volcanic hazard assessment is complex. It is recommended that at the PIA stage (possibly
even at the LS stage, depending on the site) NUMO utilize an expert elicitation process forassessment of volcanic hazards, using national and perhaps international experts.

5.3

Closing Remarks

Clearly, we have only broken into the surface of a detailed, process-based understanding of
the tectonics of Tohoku. With further development and more detailed application of the82methodology, using more extensive local expert input, we know that much more detail could
be obtained. If volunteer sites do emerge in northern Honshu, there will be an opportunity to
deploy the methodology in a more focussed and intensive manner.
Finally, and very significantly for NUMO and the national and international geosciences
community that is tracking its work, we would note that, based on the Tohoku case study, it
appears that a large number of potentially suitable sites could be found in similar, active
tectonic environments in Japan.

83

6 References
Altamimi, Z., P. Sillard, and C. Boucher (2002). ITRF2000: A new release of the International
Terrestrial Reference Frame for earth science applications, J. Geophys. Res., 107(B10), 2214,
doi:10.1029/2001JB000561.
Beanland, S., G. H. Blick, and D.J. Darby (1990). Normal faulting in a back-arc basin:
Geological and geodetic characteristics of the 1987 Edgecumbe earthquake, New Zealand, J.
Geophys. Res., 95, 4693-4707.
Beavan, J., and, J. Haines (2001). Contemporary horizontal velocity and strain-rate fields of
the Pacific-Australian plate boundary zone through New Zealand, J. Geophys. Res. 106, 741770.
Beavan, J., P. Tregoning, M. Bevis, T. Kato, and C. Meertens (2002). The motion and rigidityof the Pacific Plate and implications for plate boundary deformation, J. Geophys. Res., 107,
2261, doi:10.1029/2001JB000282.
Bennett, R., A.M. Friedrich, and K.P. Furlong (2004). Codependent histories of the San
Andreas and San Jacinto fault zones from inversion of fault displacement rates, Geology,
32(11), 961-964.
Beutler, G., H. Bock, E. Brockmann, R. Dach, P. Fridez, W. Gurtner, U. Hugentobler, D.
Ineichen, J. Johnson, M. Meindl, L. Mervart, M. Rothacher, S. Schaer, T. Springer, R. Weber(2001). Bernese GPS Software Version 4.2, Ed. by U. Hugentobler, S. Schaer, P. Fridez,
Astronomical Institute, University of Berne.
Bird, P. (2003). An updated digital model of plate boundaries, Geochem. Geophys. Geosyst.,
4(3), 1027, doi:10.1029/2001GC000252.
Boullier, A-M., Fujimoto, K., Ito, H., Ohtani, T., Keulen, N., Fabbri, O., Amitrano, D., Dubois,
M., Pezard, P. (2004). Structural evolution of the Nojima fault (Awaji Tsland, Japan) revisitedfrom the GSJ drill hole at Hirabayashi. Earth Planets Space, 56, 1233-1240.
Calais E., M. Vergnolle, V. San'kov, A. Lukhnev, A. Miroshnitchenko, S. Amarjargal, J.
Déverchère (2003). GPS measurements of crustal deformation in the Baikal-Mongolia area
(1994–2002): Implications for current kinematics of Asia, J. Geophys. Res., 108 (B10), 2501,
doi:10.1029/2002JB002373.
Chapman, N. A., Goto, J. and Tsuchi, H. (2009 in press). Likelihood of tectonic activity
affecting the geological stability of a repository in Japan: development of MUNO’s ITM
methodology. In: Stability and buffering capacity of the geosphere for long-term isolation ofradioactive waste: application to crystalline rock. Proceedings of an IGSC Geosphere Stability
Workshop, September 2007. OECD Nuclear Energy Agency, Paris.
Committee for the Catalogue of Quaternary Volcanoes in Japan (1999): Catalogue of
Quaternary Volcanoes in Japan, v 1.0 (CD-ROM), Volcanol. Soc. Japan.
Connor, C. B. and B. E. Hill (1995). Three nonhomogeneous Poisson models for the
probability of basaltic volcanism: Application to the Yucca Mountain region. Journal ofGeophysical Research, 100, 10,107−10,125.
Connor, C. B., S. Magsino and J. Stamatakos et al. (1997). Magnetic surveys help reassess
volcanic hazards at Yucca Mountain. Eos, Transactions of the American Geophysical Union,
78(7), 73, 77−78.
Connor, C. B., J. Stamatakos, D. Ferrill, B. E. Hill, G. Ofoegbu and F. M. Conway (2000).
Volcanic hazards at the proposed Yucca Mountain, Nevada, high-level radioactive waste
repository. Journal of Geophysical Research, 105, 417−432.

84Connor C., L. Connor. (2008). Estimating spatial intensity with kernel methods, In: Connor, C.,
Chapman, N. A. and Connor, L. eds. (2008, in press). Volcanic and Tectonic Hazard
Assessment for Nuclear Facilities. Cambridge University Press.
Cox D. R. (1955). Some statistical methods connected with series of events, Journal of theRoyal Statistical Society, Series B (Methodological), XVII, No. 2, pp. 129-164.
Cox, A., and Engebretson, D. (1985). Change in motion of Pacific plate at 5 Myr BP: Nature,313, 472-474.
Cox, A. and R.B. Hart (1986). Plate tectonics: How it works, Palo Alto: Blackwell Scientificpublications, 392 pp.
Crowe, B. M., D. T. Vaniman and W. J. Carr (1983). Status of volcanic hazard studies for the
Nevada nuclear waste storage investigations, Lab Report LA, LA-9325-MS. Los Alamos
National Lab, NM.
Davidson, J., and De Silva, S. (2000). Composite volcanoes: in Sigurdsson, H., ed.,Encyclopedia of Volcanoes: Academic Press, San Diego, p. 663-681.
DeMets, C., 1992, Oblique convergence and deformation along the Kuril and Japan Trenches, J.
Geophys. Res., 97(B12), 17,615-17,625.
DeMets, C., R.G. Gordon, D.F. Argus, and S. Stein (1994). Effect of recent revisions to thegeomagnetic reversal time scale on estimates of current plate motions, Geophys. Res. Lett.,
21, 2191-2194.
Diggle, P. (1985). A kernel method for smoothing point process data. Applied Statistics, 34,
138−147.
Diggle, P. and J. S. Marron (1988). Equivalence of smoothing parameter selectors in density
and intensity estimation. Journal of the American Statistical Association, 83, 793−800.
Dixon, T. H., A. Mao, M. Bursik, M. Heflin, J. Langbein, R. Stein, and F. Webb (1997).
Continuous monitoring of surface deformation at Long Valley Caldera, California, with GPS, J.
Geophys. Res., 102, 12,017– 12,034.
Dixon, T.H., E. Norabuena, and L. Hotaling (2003). Paleoseismology and Global Positioning
System: Earthquake-cycle effects and geodetic versus geological fault slip rates in the
Eastern California Shear Zone, Geology, 31(1), 55-58.
Dragert, H., K. Wang, and T. James (2001). A silent slip event on the deeper Cascadia
subduction interface, Science, 292, 1525-1528.
Engebretson, D. C., Cox, A., and Gordon, R. G. (1985). Relative motions between oceanicand continental plates in the Pacific basin: Geological Society of America Special Paper 206,
55 pp.
Finn, C. (1994). Aeromagnetic evidence for a buried early Cretaceous magmatic arc,
northeast Japan: Journal of Geophysical Research, 99, 22,165-22,185.
Finn, C., Kimura, G., Suyehiro, K. (1994). Introduction to the special section Northeast
Japan: A case history of subduction: Journal of Geophysical Research, 99, 22,137-22,145.
Frankel, A.D., Petersen, M.D., Mueller, C.S., Haller, C.M., Wheeler, R.L., Leyendecker, E.V.,
Wesson, R.L., Harmsen, S.C., Cramer, C.H., Perkins, D.M. and Rukstales, K.S. (2002).

85Documentation for the 2002 Update of the National Seismic Hazard Maps. U.S. Geological
Survey Open-File Report 02-420.
Gatrell, A. C., T. C. Bailey, P. J. Diggle and B. S. Rowlingson (1996). Spatial point patternanalysis and its application in geographical epidemiology. Transactions of the Institute ofBritish Geography NS 21, 256−274.
Haines, A. J., W. E. Holt (1993). A procedure for obtaining the complete horizontal motions
within zone of distributed deformation from the inversion of strain rate data, J. Geophys. Res.,
98(B7), 12057-12082, 10.1029/93JB00892.
Hasegawa, A., S. Horiuchi and N. Umino (1994). Seismic structure of the northeastern Japan
convergent margin: A synthesis, J. Geophys. Res., 99(B11), 22295-22312,
10.1029/93JB02797.
Hasegawa, A., A. Yamamoto, N. Umino, S. Miura, S. Horiuchi, D. Zhao, and H. Sato (2000).
Seismic activity and deformation process of the overriding plate in the northeastern Japansubduction zone, Tectonophysics, 319(4), 225-239.
Hasegawa, A., J. Nakajima (2004). Geophysical constraints on slab subduction and arc
magmatism, Geophysical monograph 150, IUGG, Volume 19, pp. 81-93.
Hasegawa, A., J. Nakajima, N. Umino, S. Miura (2005). Deep structure of the northeasternJapan arc and its implications for crustal deformation and shallow seismic activity,
Tectonophysics, 403, pp. 59-75.
Heki, K. (2004). Space geodetic evidence of deep basal subduction erosion in northeasternJapan: Earth and Planetary Science Letters, 219, 13-20.
Heki, K., S. Miyazaki and H. Tsuji (1997). Silent fault slip following an interplate thrustearthquake at the Japan Trench, Nature, 386, 595-597.
Heki, K., S. Miyazaki, H. Takahashi, M. Kasahara, F. Kimata, S. Miura, N. Vasilenco, A.
Ivashchenco and K. An (1999). The Amurian plate motion and current plate kinematics inEastern Asia, J. Geophys. Res., 104, 29147-29155.
Herring, T. A. (2001). GLOBK global Kalman filter VLBI and GPS analysis program, version
5.03, Mass. Inst. of Technol., Cambridge.
Hirose, H., K. Hirahara, F. Kimata, N. Fujii, and S. Miyazaki (1999). A slow thrust slip event
following the two 1996 Hyuganada earthquakes beneath the Bungo Channel, southwest
Japan, Geophys. Res. Lett., 26(21), 3237-3240.
Honda, S. (1985). Thermal structure beneath Tohoku, northeast Japan -- A case study forunderstanding the detailed thermal structure of the subduction zone: Tectonophysics, 112,
69-102.
Hone, D.W.E., S. Mahony, R.S.J. Sparks and K. Martin (2007). Cladistics analysis applied tothe classification of volcanoes. Bulletin of Volcanology, 70, 203-220.
Hsu, Y-J., N. Bechor, P. Segall, S.-B. Yu, L.-C. Kuo, and K.-F. Ma (2002). Rapid afterslip
following the 1999 Chi-Chi, Taiwan Earthquake, Geophys. Res. Lett., 29(16),doi:10.1029/2002GL014967.
Jaquet, O. and R. Carniel (2006). Estimation of volcanic hazards using geostatistical methods.
In: H. M. Mader et al. (eds), Statistics in Volcanology, Special Publications of IAVCEI 1.
London: Geological Society, 77−87.
Jaquet, O., C. B. Connor and L. Connor (2008). Long-term volcanic hazard assessments fornuclear facilities. Nuclear Technology, 163, 180−189.

86Jaquet O., C. Lantuéjoul (2008). Cox process models for the estimation of long-term volcanic
hazard, In: Connor, C., Chapman, N. A. and Connor, L. eds. (2008, in press). Volcanic andTectonic Hazard Assessment for Nuclear Facilities. Cambridge University Press.
Jarvis A., H. I. Reuter, A. Nelson and E. Guevara (2006). Hole-filled seamless SRTM data V3,
International Centre for Tropical Agriculture (CIAT), available from: http://srtm.csi.cgiar.org.
Jolivet, L, K., Tamaki, K., Fournier, M. (1994). Japan Sea, opening history and mechanism:
A synthesis: Journal of Geophysical Research, 99, 22,237-22,259.
Katz, Y., Weinberger, R., Aydin, A. (2004). Geometry and kinematic evolution of reidal shear
structures, Capitol Reef National Park, Utah. Journal of Structural Geology, 26, 491-501.
Kawakatsu, H., and Seno, T. (1983). Triple seismic zone and regional variation of seismicity
along the northern Honshu arc: Journal of Geophysical Research, 88, 4215-4230.
Kondo, H., Tanaka, K., Mizuochi, Y., and Ninomiya, A. (2004). Long-term changes in
distribution and chemistry of middle Miocene to Quaternary volcanism in the Chokai-Kurikoma
area across the Northeast Japan Arc: Island Arc, 13, 18-46.
Kondo H. (2009). Regional-scale volcanology in support of site-specific ienvestigations: anempirical approach based on spatio-temporal patterns in volcanism combined with the related
phenomena, In: Connor, C., Chapman, N. A. and Connor, L. eds. (2009). Volcanic and Tectonic
Hazard Assessment for Nuclear Facilities. Cambridge University Press.
Kostrov, B.V. (1974). Seismic moment and energy of earthquakes and seismic flow of rock.
Izv, Acad, Sci. USSR Physics. Solid Earth 1, 23-40.
Kreemer, C., W. E. Holt, S. Goes, R. Govers (2000). Active deformation in eastern Indonesia
and the Philippines from GPS and seismicity data, J. Geophys. Res., 105(B1), 663-680,
10.1029/1999JB900356.
Langbein J. (2004). Noise in two-color electronic distance meter measurements revisited, J.
Geophys. Res., 109, B04406, doi:10.1029/2003JB002819.
Lantuéjoul C. (2002). Geostatistical simulation: models and algorithms. Springer, 256 pp.
Larson, K., A. Lowry, V. Kostoglodov, W. Hutton, O. Sanchez, K. Hudnut, and G. Suarez(2004). Crustal deformation measurements in Guerrero, Mexico, J. Geophys. Res., 109(B4),
B04409.
Lutz, T. M. and J. T. Gutmann (1995). An improved method for determining alignments ofpoint-like features and its implications for the Pinacate volcanic field, Sonora, Mexico. Journalof Geophysical Research, 100, 17 659−17 670.
Mahoney S., R.S.J. Sparks (2008). Assessing long-term volcanic hazards using volcano
distribution databases, In: Connor, C., Chapman, N. A. and Connor, L. eds. (2008, in press).
Volcanic and Tectonic Hazard Assessment for Nuclear Facilities. Cambridge University Press.
Mansinha, L., and D.E. Smylie (1971). The displacement fields of inclined faults, Bull. Seism.
Soc. Am., 61(5), 1433-1440.
Martin, A. J., K. Umeda, C. B. Connor, J. N. Weller, D. Zhao and M. Takahashi (2004).
Modeling long-term volcanic hazards through Bayesian inference: example from the Tohoku
volcanicarc,
Japan.
JournalofGeophysicalResearch,
109,
B10208,
DOI:10.1029/2004JB003201.

87

Matsuda T. (1977). Estimation of future destructive earthquakes from active faults on land inJapan, J. Phys. Earth, Suppl. 25, 795-855.
Matsuda, T., H. Yamazaki, T. Nakata, and T. Imaizumi (1980). The surface faults associatedwith the Rikku earthquake of 1896, Bull. Earthq. Res. Inst., Univ. of Tokyo, 55, 795-855.
Mattioli, G. S., T. H. Dixon, F. Farina, E. S. Howell, P. E. Jansma, and A.L. Smith (1998).
GPS measurement of surface deformation around Soufriere Hills Volcano, Montserrat, from
October 1995 to July 1996, Geophys.Res. Lett., 25, 3417– 3420.
Mazzotti, S., X. Le Pichon, P. Henry, and S. Miyazaki (2000). full interseismic locking of the
Nankai and Japan-west Kurile subduction zones: An analysis of uniform elastic strain
accumulation in Japan constrained by permanent GPS, J. Geophys. Res., 105(B6), 13,15913,177.
McCaffrey, R. (1995). DEF-NODE users guide, Rensselaer Polytechnic Institute, Troy, New
York, www.rpi.edu/~mccafr/defnode.
McCaffrey, R. (2002). Crustal block rotations and plate coupling, in Plate Boundary Zones, pp.
100-122, eds. S. Stein and J. Freymueller, AGU Geodynamics Series, 30.
McCaffrey, R. (2005). Block kinematics of the Pacific-North America plate boundary in thesouthwestern United States from inversion of GPS, seismological and geologic data, J.
geophys. Res., 110, doi:10.1029/2004JB003307.
McCaffrey, R., M.D. Long, C. Goldfinger, P.C. Zwick, J.L. Nabelek, C.K. Johnson, and C.
Smith (2000). Rotation and plate locking at the southern Cascadia subduction zone,
Geophysical Research Letters, 27, 3117-3120.
McClusky, S. et al. (2001). Present day kinematics of the eastern California Shear zone froma geodetically constrained block model, Geophys. Res. Lett., 28, 3369-3372.
Meade, B.J., and Hager, B.H. (2005). Block models of crustal motion in southern California
constrainedbyGPS
measurements,
J.
Geophys.
Res.,
110,
B03403,
doi:10.1029/2004JB003209.
Miura, S., S. Ueki, T. Sato, K. Tachibana, and H. Hamaguchi (2000). Crustal deformation
associated with the 1998 seismo-volcanic crisis of Iwate Volcano, northeastern Japan, as
observed by a dense GPS network, Earth Planets Space, 52, 1003-1008.
Miura, S., T. Sato, A. Hasegawa, Y. Suwa, K. Tachibana and S. Yui (2004a). Strain
concentration zone along the volcanic front derived by GPS observations in the NE Japan arc,
Earth Planets Space, 56, 1347-1355.
Miura, S., Y. Suwa, A. Hasegawa, and T. Nishimura (2004b). The 2003 M8.0 Tokachi-Oki
earthquake—How much has the great event paid back slip debts?, Geophysical Res. Letters,
31, doi:10.1029/2003GL019021.
Nakata, T., Imaizumi, T. (2002). Digital active fault map of Japan.
Nikolaidis, R. (2002). Observation of geodetic and seismic deformation with the Global
Positioning System, Ph.D. thesis, Univ. of Calif., San Diego.
Nishimura, T., T. Hirasawa, S. Miyazaki, T. Sagiya, T. Tada, S. Miura, and K. Tanaka (2004).
Temporal change of interplate coupling in northeastern japan during 1995-2002 estimatedfrom continuous GPS observations, Geophys. J. Int., 157, 901-916.

88

Ohtake, M., Taira, A., Ota, Y. (2002). Active faults and seismo-tectonics of the eastern marginof the Japan Sea. University of Tokyo Press.
Okada, Y. (1985). Surface deformation due to shear and tensile faults in a half-space, Bull.
Seism. Soc. Am., 75, 1135-1154.
Okamura, Y., Watanabe, M., Morijiri, R., Satoh, M. (1995). Rifting and basin inversion in theeastern margin of the Japan Sea: Island Arc, 4, 166-181.
Okubo, Y., and Matsunaga, T. (1994). Curie-point depth in northeast Japan and its
correlation with regional thermal structure and seismicity: Journal of Geophysical Research,
99, 22,363-22,371.
Ota, Y., Omura, A. (1991). Late Quaternary shorelines in the Japanese Islands. DaiyonkiKenkyu=The Quaternary Research, 30,175–186.
Parsons, T., G. A. Thompson and A. H. Coghill (2006). Earthquake and volcano clustering via
stress transfer at Yucca Mountain, Nevada. Geology, 34, 785−788, doi:10.1130/G22636.1.
Pollitz, F. F. (1986). Pliocene change in Pacific-plate motion: Nature, 320, 738-741.
Pollitz, F. (1997). Gravitational viscoelastic postseismic relaxation on a layered spherical
Earth, J. Geophys. Res., 102(B8), 17,921-17,941.
Pollitz, F., M. Nyst, T. Nishimura, and W. Thatcher (2006). Inference of postseismic
deformation mechanisms of the 1926 Kanto earthquake, J. Geophys. Res., 111,
doi:10.1029/2005JB003901.
Prawirodirdjo, L. and Y. Bock (2004). Instantaneous global plate motion model from 12 yearsof continuous GPS observations, J. Geophys. Res., 109, B08405, doi:10.1029/2003JB002944.
Pyle, D. M. (2000). Sizes of volcanic eruptions. In: H. Sigurdsson et al. (eds), Encyclopedia ofVolcanoes. New York: Academic Press, 263−271.
Reilinger R., et al. (2006). GPS constraints on continental deformation in the Africa-ArabiaEurasia continental collision zone and implications for the dynamics of plate interactions, J.
Geophys. Res., 111, B05411, doi:10.1029/2005JB004051.
Reiter, L., (1990). Earthquake hazard analysis. Columbia University Press. 254pp.
Rothacher, M., and L. Mervart (eds.) (1996). Documentation of the Bernese GPS Software
Version 4.0, 418 pp, Astron. Inst., Univ. of Bern, Bern, Switzerland.
Sagiya, T. and W. Thatcher (1999). Coseismic slip resolution along a plate boundary
megathrust: The Nankai Trough, southwest Japan, J. Geophys. Res., 104(B1), 1111-1129.
Sagiya, T., S. Miyazaki, and T. Tada (2000). Continuous GPS array and present-day crustal
deformation of Japan. Pure and Applied Geophysics, 2303-2322.
Sato, H. (1994). The relationship between Late Cenozoic tectonic events and stress field
and basin development in northeast Japan: Journal of Geophysical Research, 99, 22,26122,274.
Sato, H., and Amano, K. (1991). Relationship between tectonics, volcanism, sedimentation
and basin development, late Cenozoic, central part of northern Honshu, Japan: Sedimentary
Geology, 74, p. 323-343.

89

Savage, J. C., J. L. Svarc (1997). Postseismic deformation associated with the 1992 Mw=7.3
Landers earthquake, southern California, J. Geophys. Res., 102(B4), 7565-7578,
10.1029/97JB00210.
Segall, P., R. Bürgmann, M. Matthews (2000). Time-dependent triggered afterslip following
the 1989 Loma Prieta earthquake, J. Geophys. Res., 105(B3), 5615-5634,
10.1029/1999JB900352.
Sella, G.F., T.H. Dixon, and A.L. Mao (2002). REVEL: A model for recent plate velocities fromspace geodesy: J. Geophys. Res., 107(B4), doi:10.1029/2000JB000033.
Senior Seismic Hazard Analysis Committee (SSHAC) (1997). Recommendations forprobabilistic seismic hazard analysis: Guidance on uncertainty and use of experts.
NUREG/CR-6372. UCRL-ID-122160
Seno, T., and T. Sakurai, S. Stein (1996). Can the Okhotsk plate be discriminated from theNorth American plate, J. Geophys. Res., 101(B5), 11,305-11,315.
Sibson, R.H. (2003). Thickness of the seismic slip zone. Bulletin of the Seismological Societyof America, 93, 1169-1178.
Stepp, J.C., Wong, I., Whitney, J., Quittmeyer, R., Abrahamson, N., Toro, G., Youngs, R.,
Coppersmith, K., Savy, J., Sullivan, T., and Yucca Mountain PSHA project members (2001).
Probabilistic seismic hazard analyses for ground motions and fault displacement at Yucca
Mountain, Nevada. Earthquake Spectra, 17, 113-151
Stirling, M.W., Wesnousky, S.G., Shimazaki, K. (1996). Fault trace complexity, cumulative slip,
and the shape of the magnitude-frequency distribution for strike-slip faults: a global survey
Geophysical Journal International, 124 (3) , 833–868.
Stirling, M.W., McVerry, G.H., and Berryman, K.R. (2002). A new seismic hazard model for
New Zealand, Bull. Seism. Soc. Am., 92, 1878-1903.
Subarya, C., Chlieh, M., Prawirodirdjo L., Avouac J-P., Bock, Y., Sieh, K., Meltzner, A.J.,
Natawidjaja, D.H., McCaffrey, R. (2006). Plate-boundary deformation associated with thegreat Sumatra–Andaman earthquake. Nature, 440, 46-51.
Sugimura, A., and S. Uyeda (1973). Island Arcs: Japan and its environs, Elsevier, Amsterdam,Suwa, Y., S. Miura, A. Hasegawa, T. Sato, and K. Tachibana (2006). Interplate coupling
beneath NE Japan inferred from three-dimensional displacement field, 2006, J. Geophys.
Res., 111, doi:10.1029/2004JB003203.
Swofford, D. L. (1998). PAUP* Phylogenetic analysis using parsimony (*and other methods),
version 4. Sinauer, Massachusetts.
Sylvester, A.G. (1988). Strike-slip faults. Geological Society of America Bulletin, 100, 16661703.
Taira, A. (2001). Tectonic evolution of the Japanese island arc system: Annual Reviews ofEarth and Planetary Sciences, 29, p. 109-134.
Tajikara, M., Ikeda, Y. (2005). Vertical Crustal Movement and Development of Basin andRange Topography in the Middle Part of the Northeast Japan Arc Estimated from
Fluvial/Marine Terrace Data. Quaternary Research, 44, 229-245.
Takahashi, H. and 13 other authors (1999). Velocity field of around the Sea of Okhotsk andSea of Japan regions determined from a new continuous GPS network data, Geophys. Res.
Lett., 26(16), 2533-2536.

90

Takeuchi, H., Uyeda, S., and Kanamori, H. (1970). Debate about the Earth: Approach toGeophysics through Analysis of Continental Drift: Freeman, Cooper and Company, San
Francisco, 281 pp.
Tamaki, K., and Honza, E. (1985). Incipient subduction and obduction along the eastern
margin of the Japan Sea: Tectonophysics, 119, p. 381-406.
Tamura Y., Y. Tatsumi, D. Zhao, Y. Kido, H. Shukuno (2002). Hot fingers in the mantle
wedge: new insights into magma genesis in subduction zones, Earth and Planetary Sciences
Letters, 197, pp. 105-116.
Tamura Y., J. Nakajima, S. Kodaira (2008). Tectonic setting of volcanic centers in subduction
zones: 3-D structure of mantle wedge and arc crust, In: Connor, C., Chapman, N. A. andConnor, L. eds. (2008, in press). Volcanic and Tectonic Hazard Assessment for Nuclear
Facilities. Cambridge University Press.
Tchalenko, J.S., Ambraseys, N.N. (1970). Structural analysis of the Dasht-e Bayaz (Iran)earthquake fractures. Geological Society of America Bulletin, 81, 41-60.
Thatcher, W. and J. B. Rundle (1984). A viscoelastic coupling model for the cyclic
deformation due to periodically repeated earthquakes at subduction zones, J. Geophys. Res.,89, 7631-7640.
Von Huene, R., and Lallemand, S. (1990). Tectonic erosion along the Japan and Peru
convergent margins: Geological Society of America Bulletin, 85, p. 1159-1170.
Wald, D.J., and Graves, R.W. (2001). Resolution analysis of finite fault source inversion
using one- and three-dimensional Green’s functions. 2. Combining seismic and geodetic
data. J. Geophysical Research, 106, 8767-8788.
Wallace, K., G. Yin, and R. Bilham (2004). Inescapable slow slip on the Altyn Tagh Fault,
Geophysical Research Letters, 31, doi:10.1029/2004GL019724
Wallace, L.M., J. Beavan, R. McCaffrey and D. Darby (2004). Subduction zone coupling andtectonic block rotations in the North Island, New Zealand, J. Geophys. Res., 109(B12), doi:
10.1029/2004JB003241.
Wallace, L.M., J. Beavan, R. McCaffrey, K. Berryman, and P. Denys (2007). Balancing theplate motion budget in the South Island, New Zealand using GPS, geological andseismological data, Geophys. J. Int., doi: 10.1111/j.1365-246X.2006.03183.x
Wei, D. and T. Seno (1998). Determination of the Amurian Plate motion, in, Mantle Dynamics
and Plate Interactions in East Asia, AGU Geodynamics, 27, 337-346.
Weichert, D.H. (1980). Estimation of the earthquake recurrence parameters for unequalobservation periods for different magnitudes. Bulletin of the Seismological Society of America,
70, 1337-1346.
Wesnousky, S.G., Scholz, C.H., Shimazaki, K., Matsuda, T. (1984). Integration of geologicaland seismological data for the analysis of seismic hazard: A case study of Japan. Bulletin of
the Seismological Society of America, 74, 687-708.
Wetmore, P. H., S. S. Hughes, L. J. Connor, and M. L. Caplinger (2008 in press). Spatial
distribution of eruptive centers about the Idaho National Laboratory. In: Connor, C. B.,
Chapman, N. A., Connor, L. J. (eds) Volcanic and Tectonic Hazard Assessment for Nuclear
Facilities. Cambridge, UK: Cambridge University Press.
Williams S. D. P., Y. Bock, P. Fang, P. Jamason, R. M. Nikolaidis, L. Prawirodirdjo, M. Miller,
D. J. Johnson (2004). Error analysis of continuous GPS position time series, J. Geophys.
Res., 109, B03412, doi:10.1029/2003JB002741.

91

Woo G. (1999). The mathematics of natural catastrophes. Imperial College Press, London,
292 pp.
Yamanaka Y., and M. Kikuchi (2004). Asperity map along the subduction zone innortheastern Japan inferred from regional seismic data, J. Geophys. Res., 109, B07307,
doi:10.1029/2003JB002683.
Zhang, J., Y. Bock, H. Johnson, P. Fang, S. Williams, J. Genrich, S. Wdowinski, and J. Behr
(1997). Southern California Permanent GPS Geodetic Array: Error analysis of daily position
estimates and site velocities, J. Geophys. Res., 102(B8), 18,035–18,056.
Zhao, D. (2001). Seismological structure of subduction zones and its implications for arc
magmatism and dynamics. Physics of the Earth and Planetary Interiors, 127, 197−214.
Zhao D., A. Hasegawa, S. Horiuchi (1992). Tomographic imaging of P and S wave velocity
structure beneath Northeastern Japan, Journal of Geophysical Research, 97, pp. 1990919928.
Zhao D., A. Hasegawa, H. Kanamori (1994). Deep structure of Japan subduction zone asderived from local, regional and teleseismic events, Journal of Geophysical Research, 99, pp.
22313-22329.
Zhao D., F. Ochi, A. Hasegawa, A. Yamamoto (2000). Evidence for the location and cause oflarge crustal earthquakes in Japan, Journal of Geophysical Research, 105, pp. 13579-13594.
Zhu L. and Rivera, L.A. (2002). A note on the dynamic and static displacements from a point
source in multilayered media. Geophys. J. Int. 148, 619-627.

92

Appendix
Descriptions of the Fourteen Example Locations
in Tohoku and the Results of the Probabilistic
Evaluations for Each Location

Background
The overall tectonic setting of northern Honshu and Tohoku is described in Section 3 of thisreport.
Northern Honshu is within the plate boundary zone between the Pacific and Amurian plates,
which are converging at ~9-10 cm/yr (e.g., DeMets, 1992; Seno et al., 1996; Heki et al., 1999).
Westward subduction of the Pacific Plate beneath northern Honshu occurs at the Japan
Trench. The active volcanic arc related to this subduction occurs within central northern
Honshu. The location of the arc in this region of Japan has been approximately stable since
the late Pliocene (see review in Finn et al., 1994). Although most (~80%) of the relative platemotion is accommodated on the Japan Trench subduction thrust (e.g., DeMets, 1992; Seno etal., 1996), an active convergent zone in the eastern Sea of Japan accommodates
underthrusting of the Sea of Japan crust beneath northern Honshu (up to 1-2 cm/yr of
convergence; e.g., Seno et al., 1996; Heki et al., 1999) via a series of folds and thrusts
offshore, and in western northern Honshu. It is suggested that this contractional zone may
evolve into a subduction zone in the future (Tamaki and Honza, 1985). Smaller amounts ofcontractional deformation also occur in the Ou/Backbone range (e.g., Miura et al., 2004) ofcentral northern Honshu via high-angle reverse faults. The Ou/Backbone Range zone ofcontraction coincides with the main volcanic arc (Hasegawa et al., 2000). Between theeastern margin of the Ou/Backbone range and the east coast of the Tohoku region is the
Kitakami Range, a forearc basement block that is considered to have been tectonically stablesince the Miocene (Sato and Amano, 1991), although it has undergone some uplift in the lateCenozoic (Sugimura and Uyeda, 1973).
(Citations are included in the main reference list of the report).

Example Location 1: 39.480N 140.490E
Tectonic setting
Location 1 is situated between the Backbone Range and the fold and thrust zone related to
the Sea of Japan convergence zone. The example location is on the western side of the
footwall basin of the Senya Fault, which is located 12 km to the east. The uplift contours show
long-term tilting toward the fault (east). The example location is 50 km south-west of Sengan,
the nearest volcano is Inai-Takahachi at 35 km to the north-east and lies within a volcanic‘gap’ region.

Geological environment
Location 1 is situated within Quaternary sand, gravel and clay of the footwall basin. The total
thickness example location is unknown, but thins westward, so may be only tens of metres.
The bedrock beneath the basin is likely to be late Miocene mudstone (with sandstone andacid tuff) and middle Miocene sandstone, mudstone, acid pyroclastic rock and conglomerate,
which outcrops to the west.

Individual strain indicator results
Derivation of the discrete probability histograms shown below and cumulative probability
histogram (next section), and the future time periods over which the strain rates are
considered valid are described in Section 3.

-1-

Location 1 - GPS
0.4
0.35

Weighted Frequency

0.3
0.25
0.2
0.15
0.1
0.05
0
1

5

9 13 17 21 25 29 33 37 41 45 49 53 57 61 65 69 73 77 81 85 89
nanotrains/year

Location 1 - SEIS
0.4
Weighted Frequency

0.35
0.3
0.25
0.2
0.15
0.1
0.05
1
5
9
13
17
21
25
29
33
37
41
45
49
53
57
61
65
69
73
77
81
85
89

0
nanos trains /ye ar

Location 1 - SURF DEF
0.4

Weighted Frequency

0.35
0.3
0.25
0.2
0.15
0.1
0.05
0
1 5

9 13 17 21 25 29 33 37 41 45 49 53 57 61 65 69 73 77 81 85 89

nanos trains /yr-2-

Combined strain budgets and probabilities

Location 1

1
Probability of Exceedance

GPS
SEIS

0.8

SURF DEF
EQ W GTD AV

0.6
0.4
0.2
0
1

10

100

1000

nanostrains/year
GPS: High probability of exceedance of 10 nanostrain/yr (60%), but low probability ofexceedance of 20 nanostrain/yr (20%).
SEIS: High probability of exceedance of 5 nanostrain/yr (60%), but low probability ofexceedance of 8 nanostrain/yr (20%).
SURFACE DEFORMATION: High probability of exceedance of 25 nanostrain/yr (60%), but
low probability of exceedance of 27 nanostrain/yr (20%).
STRAIN RATE UNCERTAINTY
Percentile

GPS

Surf Def

Seis

Weight

16th

2

22

4

3

50th

3

24

10

11

84th

7

26

18

24

The seismicity dataset is incomplete because the earthquake sequence associated with the1896 Rikyuu Earthquake is missing from catalogue. Since the example location is within 12
km of an active fault, which has a long recurrence interval, the seismicity dataset may
underestimate tectonic strain in the area. The surface deformation dataset is likely to becomplete because the example location is situated in an area with river terraces. The GPS
dataset may be reflecting strain related to interseismic coupling processes on nearby faults(e.g., within 10-20 km) or may represent strain due to faulting very close to, or at the example
location. The dataset is consistent with the example location being within a contractionalstrain regime, and this strain may be reflecting the presence of hidden reverse faults near theexample location. The average curve gives a high probability of exceedance of 9
nanostrain/yr (60%), but low probability of exceedance of 25 nanostrain/yr (20%).

-3-

Magma intrusion probabilities and volcanic/tectonic strain interactions
Volcano-tectonic interactions at this location are unlikely, because the nearest volcanic edifice
is 30-50 km distant, and other phenomena that have been linked to volcanism, such
underplating of magmatic material leading to elevated plateau are lacking. This location is 55
km south-west of the centre of the Sengan cluster of volcanoes, the nearest volcano is InaiTakahachi at 35 km to the north-east. It is in the lowlands, between the higher ground of the
volcanic front to the east and the area of back-arc volcanism to the west, approximately half
way between clusters 5 and 6, with no known large active faults nearby.
Locations that are 30-70 km from the Sengan cluster all show a similar median, range andmaximum, log-mean and arithmetic mean with no tendency for probabilities to change with
distance. The hazard is about 14 to 30% of the hazard at location 3. The alternate dataset
analysis for this particular location (using dataset 3 as a reference dataset) suggests it has
fairly low probabilities for new volcanism in the next 1 Ma, with the logarithmic mean,arithmetic mean and maximum probabilities being 0.016, 0.0067 and 0.089 respectively. Thislocation was chosen because it is in a low-lying area, no-where near any faults or volcanoes.

Site investigation needs and siting confidence assessment
Detailed example location geology investigations (mapping, geomorphology, structural
analysis of bedrock). Geophysical techniques (GPR, Seismics, Gravity). Determination ofrates and mechanisms of uplift using geochronology, geomorphology, numerical modelling.
Refine tilt and folding deformation by improved Quaternary mapping and dating of landforms.
If potentially active faults are found then more detailed paleoseismic techniques (trenching,drilling, surveying). Enhanced GPS network to decrease uncertainties in GPS strain rates.
Deployment of microseismic network to obtain detailed distribution of local seismicity
including depths, lineations, rate parameters, and focal mechanisms.
For example location 1, the high surface deformation strain rate requires investigation of
possibility of reverse faults on the west side of the basin. Detailed mapping of the Tertiary
rocks would confirm whether or not there has been a long-term record of brittle deformation in
the example location area. Due to the moderate levels of tectonic strain, this example location
could present moderate susceptibility to tectonic processes over the next 10,000 years. This
example location would be categorised as Level 2: Medium Confidence, requiring detailed
assessment of impact scenarios and a quantitative assessment of risks by the PA team.

-4-

Example Location 2: 39.550N 140.650E
Tectonic setting and geological environment
Location 2 is situated within the backbone range, on the upthrown side of the Senya Fault.
The location is 1 km east of the mapped fault trace. Bedrock faults also exist within a few
hundred metres of the location. The nearest volcanic edifice is Inai-Takahachi, at 20 km
distant and the location is within a volcanic ‘gap’ region and near to an active fault. It issituated within bedrock comprising middle Miocene pyroclastics, mudstone and sandstone
and possibly andesite (with basalt).
The individual strain indicators are shown at the end of the example description. Note thesurface deformation histogram is truncated (max strain = 336 nanostrain/yr).

Combined strain budgets and probabilities

Location 2

1

Probability of ExceedanceGPS
SURF DEF

0.8

SEIS
EQ W GTD AV

0.6
0.4
0.2
0
1

10

100

1000

nanostrains/year

GPS: High probability of exceedance of 20 nanostrain/yr (60%), but low probability ofexceedance of 40 nanostrain/yr (20%).
SEIS: High probability of exceedance of 4 nanostrain/yr (60%), but low probability ofexceedance of 10 nanostrain/yr (20%).
SURFACE DEFORMATION: High probability of exceedance of 70 nanostrain/yr (>60%), but
low probability of exceedance of 80 nanostrain/yr (20%).
STRAIN RATE UNCERTAINTY
Percentile

GPS

Surf Def

Seis

Weight

16th

9

61

3

4

50th

23

66

5

24

84th

34

80

11

66

-5-

The equal weighting curve of all three datasets provides a measure by which to compare
strain rate models among sites. The seismicity dataset is incomplete because the earthquake
sequence associated with the 1896 Rikyuu Earthquake is missing from catalogue. Since thelocation is within 5 km of an active fault, which has a long recurrence interval, the seismicity
dataset may underestimate tectonic strain in the area. The surface deformation dataset islikely to be complete because it is dominated by active fault data from the well-studied Senya
Fault. The GPS dataset may be reflecting strain related to interseismic coupling processes on
nearby faults (e.g., within 10-20 km) or may represent strain due to faulting very close to, or at
the location. The dataset is consistent with the location being within a contractional strain
regime. The weighted average curve gives a high probability of exceedance of 12
nanostrain/yr (60%), but low probability of exceedance of 70 nanostrain/yr (20%).

Magma intrusion probabilities and volcanic/tectonic strain interactions
Volcano-tectonic interactions at this location are unlikely, because the nearest volcanic edificeis Inai-Takahachi at 20 km distant, and other phenomena that have been linked to volcanism,
such underplating of magmatic material leading to elevated plateau are lacking. Location 2 isin “volcanic gap” between the clusters on the back-arc side of the volcanic front, 41 km fromthe centre of Sengan cluster. It is distal from the main volcano clusters and is in the back arc
region behind the volcanic front, positioned on the boundary between the high ground of the
volcanic front and the low-lying area directly to the west (behind) the front. Volcano formation
at location 2 is expected to be more likely than at a location in the low land on the fore-arc
side such as location 6. However, although localities 2 and 6 are approximately the same
distance from the volcanic front, location 2 has the slightly lower hazard values.
Locations that are 30-70 km from the Sengan cluster all show a similar median, range andmaximum, log-mean and arithmetic mean with no tendency for probabilities to change with
distance. The hazard is about 14 to 30% of the hazard at location 3. The alternate dataset
analysis for this particular location (using dataset 3 as a reference dataset), suggests it has
intermediate probabilities for new volcanism in the next 1 Ma, with the logarithmic mean,arithmetic mean and maximum probabilities being 0.022, 0.014 and 0.105 respectively,
making it a similar risk to Location 13. However Location 2 is located very close to a fault, so
future site stability would not be guaranteed.

Site investigation needs and siting confidence assessment
Detailed site geology investigations will be required (mapping, geomorphology, structural
analysis of bedrock). Geophysical techniques suggested include (GPR, Seismics, Gravity).
Determination of rates and mechanisms of uplift using geochronology, geomorphology,
numerical modelling are also required. Refine tilt and folding deformation by improved
Quaternary mapping and dating of landforms. If potentially active faults are found then moredetailed palaeoseismic techniques (trenching, drilling, surveying) will be needed. Enhanced
GPS network suggested to decrease uncertainties in GPS strain rates. Deployment of
microseismic network also proposed to obtain detailed distribution of local seismicity including
depths, lineations, rate parameters and focal mechanisms
For location 2, being within 1 km of the trace of the fault would require careful examination of
the accuracy of the fault mapping, complexity of the faulting and the width of the process zone.
Additional concerns are how the fault evolves over the next 10,000 years and the complexityof secondary faulting and fracturing in the hanging wall of a reverse fault.
If further investigations were undertaken, bedrock faults should be examined. The sandstone
and mudstone bedrock may also provide evidence of long-term brittle deformation. Due to theproximity to an active fault, this location could present high probabilities of some types oftectonic impact over the next 10,000 years, and would be categorised as Level 3: Limited
Confidence.

-6-

Location 2 - GPS
0.4
0.35

Weighted Frequency

0.3
0.25
0.2
0.15
0.1
0.05
0
1

5

9 13 17 21 25 29 33 37 41 45 49 53 57 61 65 69 73 77 81 85 89
nanotrains/ye ar

Location 2 - SEIS
0.4
Weighted frequency

0.35
0.3
0.25
0.2
0.15
0.1
0.05
1
5
9
13
17
21
25
29
33
37
41
45
49
53
57
61
65
69
73
77
81
85
89

0
nanos trains /ye ar

Location 2 - SURF DEF
0.4

Weighted Frequency

0.35
0.3
0.25
0.2

NB. truncated

0.15
0.1
0.05
0
1 5

9 13 17 21 25 29 33 37 41 45 49 53 57 61 65 69 73 77 81 85 89

nanos trains /yr-7-

Example Location 3: 39.770N 140.760E
Tectonic setting and geological environment
Example location 3 is within the backbone range, within the Sengan volcanic cluster. The
nearest active fault is 4 km to the southwest and the nearest volcanic centre is 4 km, so thislocation falls within a volcano exclusion zone. However, for the purposes for this exercise, wewill report and evaluate the location. It is very near the centre of the Sengan cluster, with theyoung Akita-Komagatake being the closest volcano, just 4 km to the south-east. The location
is situated within Pleistocene mudflow and alluvial deposits. The bedrock is Pleistocene
andesite lavas and pyroclastic rocks.
The individual strain indicators are shown at the end of the example description.

Combined strain budgets and probabilities

Location 3

Probability of Exceedance

1
GPS

0.8

SURF DEF
SEIS

0.6

EQ W GTD AV
0.4
0.2
0
1

10

nanostrains/year

100

1000

GPS: High probability of exceedance of 20 nanostrain/yr (60%), but low probability ofexceedance of 40 nanostrain/yr (20%).
SEIS: High probability of exceedance of 4 nanostrain/yr (60%), but low probability ofexceedance of 10 nanostrain/yr (20%).
SURFACE DEFORMATION: High probability of exceedance of 7 nanostrain/yr (60%), but low
probability of exceedance of 15 nanostrain/yr (20%).

STRAIN RATE UNCERTAINTY
Percentile

GPS

Surf Def

Seis

Weight

16th

12

4

3

4

50th

19

11

5

11

84th

39

26

13

26

-8-

The equal weighting curve of all three methods provides a measure by which to compare
strain rate models among locations. A small proportion of the seismicity dataset may be due
to volcanogenic processes, however all focal mechanisms in the small to moderatemagnitude range show reverse motion which is consistent with tectonic processes. The
surface deformation dataset is likely to be complete because it is dominated by wellconstrained active fault data. The GPS dataset may be reflecting strain related to interseismic
coupling processes on nearby faults or may be representing strain occurring at the location.
The dataset is consistent with the location within a contractional strain region. The weighted
average curve gives a high probability of exceedance of 8 nanostrain/yr (60%), but low
probability of exceedance of 20 nanostrain/yr (20%).

Magma intrusion probabilities and volcanic/tectonic strain interactions
Near-field (< 10 km) interactions at earthquake or eruptive crisis scale are likely at thislocation. Triggered seismicity in the region and on known active faults 4 km distant are highly
likely in the event of major volcanic activity. Conversely, a major earthquake could trigger anearby eruption if a magma chamber were susceptible. Over time, the interaction of
volcanism and tectonism at this location could promote the migration of fault activity or the
locus of volcanism toward the location. Location 3 is 15 km from the centre of Sengan clusteron the volcanic front and has the highest hazard for nearly all of the alternative volcano
datasets. This relatively high hazard is expected due to its location within the Sengan clusterof volcanoes, with the young Akita-Komagatake being the closest volcano, just 4 km to thesouth-east.
The test locations < 30 km from Sengan cluster show the highest hazard, decreasing with
distance from away from the volcano cluster. The inter-quartile range, median, maximum,
logarithmic mean and arithmetic mean are all highest for these localities. As expected, itshows the highest probability for nearly all of the alternative volcano datasets, simply due toits very close proximity to the Sengan cluster of volcanoes. The alternate dataset analysis forthis particular location (using dataset 3 as a reference dataset), suggests it has very highprobabilities for new volcanism in the next 1 Ma, with the logarithmic mean, arithmetic mean
and maximum probabilities being 0.070, 0.048 and 0.734 respectively. This location waschosen to find a possible maximum value of probability, to allow comparison with otherlocations.

Site investigation needs and siting confidence assessment
This location would be excluded from consideration for a subsequent PIA because it issituated within a volcanic cluster. The following investigations would be required: detailed site
geology investigations (mapping, geomorphology, structural analysis of bedrock); geophysicaltechniques (GPR, seismics, gravity); determination of rates and mechanisms of uplift usinggeochronology, geomorphology, numerical modelling; refined tilt and folding deformation
evaluations by improved Quaternary mapping and dating of landforms; if potentially active
faults are found then more detailed palaeoseismic techniques (trenching, drilling, surveying)
are also needed; enhanced GPS network to decrease uncertainties in GPS strain rates;
deployment of microseismic network to obtain detailed distribution of local seismicity including
depths, lineations, rate parameters and focal mechanisms.
For location 3, being within 4 km of the trace of the fault would require careful examination of
the accuracy of the fault mapping, complexity of the faulting and the width of the process zone.
An additional concern is how the fault evolves over the next 10,000 years. Due to theproximity to an active fault, this location could present high probabilities of some types oftectonic impact over the next 10,000 years, and would be categorised as Level 3: Limited
Confidence, if not excluded on the basis of the volcanic criterion.

-9-

Location 3 - GPS
0.4
0.35

Weighted Frequency

0.3
0.25
0.2
0.15
0.1
0.05

89

85

81

77

73

69

65

61

57

53

49

45

41

37

33

29

25

21

17

9

13

5

1

0
nanotrains/year

Location 3 - SEIS
0.4
Weighted Frequency

0.35
0.3
0.25
0.2
0.15
0.1
0.05
1
5
9
13
17
21
25
29
33
37
41
45
49
53
57
61
65
69
73
77
81
85
89

0
nanos trains /ye ar

Location 3 - SURF DEF
0.4

Weighted Frequency

0.35
0.3
0.25
0.2
0.15
0.1
0.05
0
1 5

9 13 17 21 25 29 33 37 41 45 49 53 57 61 65 69 73 77 81 85 89

nanos trains /yr- 10 -Example Location 4: 39.820N 140.350E
Tectonic setting and geological environment
Example location 4 is midway between the backbone range and the Sea of Japan
contractional zone. It is situated on an elevated plateau, linking the Sengan cluster with theJapan Sea Coast. The nearest active fault is 15 km to the west and inactive bedrock faults
have been mapped within 2-3 km. This location is outside the main reach of the Sengan
volcanoes, being 15 km from the western edge of the cluster. The nearest active volcano is
Daibutsu Volcano, 15 km distant. The location is within bedrock that comprises Cretaceous
granodiorite and early Tertiary andesite lavas and breccias.
The individual strain indicators are shown at the end of the example description.

Combined strain budgets and probabilities

Location 4
1

Probability of ExceedanceGPS
SURF DEF

0.8

SEIS
EQ W GTD AV

0.6
0.4
0.2
0
1

10

nanostrains/year

100

1000

GPS: High probability of exceedance of 4 nanostrain/yr (60%), but low probability ofexceedance of 6 nanostrain/yr (20%).
SEIS: High probability of exceedance of 3 nanostrain/yr (60%), but low probability ofexceedance of 8 nanostrain/yr (20%).
SURFACE DEFORMATION: High probability of exceedance of 1.5 nanostrain/yr (60%), but
low probability of exceedance of 1.8 nanostrain/yr (20%).

STRAIN RATE UNCERTAINTY
Percentile

GPS

Surf Def

Seis

Weight

16th

2

2

2

2

50th

5

2

3

3

84th

6

2

8

6

- 11 -

The equal weighting curve of all three methods provides a measure by which to compare
strain rate models among locations. Seismicity rates are likely to be a good representation of
long-term background seismicity rates because are away from active faults where preinstrumental earthquake sequences would not be recorded. The surface deformation dataset
is not very discriminating of strain because of it is dominated by the tilt data, which is poorly
constrained on the relatively flat and raised plateau at this location. The nearest active fault is
15 km to the west. The GPS dataset is likely to be a good representation of the regional strain.
The equal weighted average curve gives a high probability of exceedance of 1 nanostrain/yr(60%), but low probability of exceedance of 5 nanostrain/yr (20%).

Magma intrusion probabilities and volcanic/tectonic strain interactions
This example location is approximately 15 km from the nearest active fault and nearest active
volcano (Daibutsu Volcano) on the high plateau linking the Sengan volcanic cluster with theJapan Sea coast. Short-term, earthquake or eruption scale interactions are less likely to affectthe location than at location 3. Long-term interactions between deep-seated tectonic
processes, such as mantle flow, and future volcanism are likely, but would require anunderstanding of which processes are driving uplift of the plateau and its evolution. Location 4
is situated outside the main reach of the Sengan volcanoes, being 15 km from the westernedge of the cluster. This allows comparison to location 6 to see if the risks vary either side of
the Sengan cluster.
Locations that are 30-70 km from the Sengan cluster all show a similar median, range andmaximum, log-mean and arithmetic mean with no tendency for probabilities to change with
distance. The hazard is about 14 to 30% of the hazard at location 3. The alternate dataset
analysis for this particular location (using dataset 3 as a reference dataset), suggests it has
an intermediate probability for new volcanism in the next 1 Ma, with the logarithmic mean,arithmetic mean and maximum probabilities being 0.019, 0.0097 and 0.084 respectively.

Site investigation needs and siting confidence assessment
This location may be excluded from consideration for a subsequent PIA because it is justwithin the 15 km volcano exclusion zone. However, if further investigations were undertaken,
inactive (bedrock) faults need characterisation and the processes causing plateau uplift need
to be understood. The following investigations would be required: detailed site geology
investigations (mapping, geomorphology, structural analysis of bedrock); geophysicaltechniques (GPR, seismics, gravity); determination of rates and mechanisms of uplift usinggeochronology, geomorphology, numerical modelling; refinement of tilt and folding
deformation estimates by improved Quaternary mapping and dating of landforms; if potentially
active faults are found then more detailed palaeoseismic techniques (trenching, drilling,surveying) are needed; enhanced GPS network to decrease uncertainties in GPS strain rates;
deployment of microseismic network to obtain detailed distribution of local seismicity including
depths, lineations, rate parameters and focal mechanisms.
Due to the moderate levels of tectonic strain, this location could present moderate
susceptibility to tectonic processes over the next 10,000 years. It would be categorised as
Level 2: Medium Confidence, requiring detailed assessment of impact scenarios and aquantitative assessment of risks by the PA team.

- 12 -

Location 4 - GPS
0.4
0.35

Weighted Frequency

0.3
0.25
0.2
0.15
0.1
0.05

89

85

81

77

73

69

65

61

57

53

49

45

41

37

33

29

25

21

17

9

13

5

1

0
nanotrains/ye ar

Location 4 - SEIS
0.4
Weighted Frequency

0.35
0.3
0.25
0.2
0.15
0.1
0.05
1
5
9
13
17
21
25
29
33
37
41
45
49
53
57
61
65
69
73
77
81
85
89

0
nanos trains /ye ar

Location 4 - SURF DEF
0.4

NB. truncatedWeighted Frequency

0.35
0.3
0.25
0.2
0.15
0.1
0.05
0
1 5

9 13 17 21 25 29 33 37 41 45 49 53 57 61 65 69 73 77 81 85 89

nanos trains /yr- 13 -Example Location 5: 39.850N 140.610E
Tectonic setting and geological environment
Location 5 is situated on the elevated plateau 15 km west of the Sengan cluster, 10 km westof the edge of the backbone range. The nearest active fault is 15 km to SE and inactive(bedrock) faults have been mapped within 3 km. The nearest volcanic edifice is 8 km to the
WSW. The location is approximately 9 km from the nearest active volcano (Daibutsu Volcano)and is 15 km from the projected centre of the Sengan volcano cluster. Location 5 is situatedon bedrock, which comprises middle Miocene sandstone, mudstone, and acid pyroclastic
rocks.
The individual strain indicators are shown at the end of the example description.

Combined strain budgets and probabilities

Location 5
1

Probability of ExceedanceGPS
SURF DEF

0.8

SEIS
EQ W GTD AV

0.6
0.4
0.2
0
1

10

nanostrains/year

100

1000

GPS: High probability of exceedance of 6 nanostrain/yr (60%), but low probability ofexceedance of 7 nanostrain/yr (20%).
SEIS: High probability of exceedance of 4 nanostrain/yr (60%), but low probability ofexceedance of 8 nanostrain/yr (20%).
SURFACE DEFORMATION: High probability of exceedance of 7 nanostrain/yr (60%), but low
probability of exceedance of 8 nanostrain/yr (20%).

STRAIN RATE UNCERTAINTY
Percentile

GPS

Surf Def

Seis

WGT

16th

2

7

3

3

50th

6

7

4

6

84th

7

8

10

8

- 14 -

The equal weighting curve of all three methods provides a measure by which to compare
strain rate models among locations. Seismicity rates are likely to be a good representation of
long-term background seismicity rates because they are away from active faults where preinstrumental earthquake sequences would not be recorded. The surface deformation dataset
is not very discriminating of strain because of it is dominated by the tilt data (the nearest
active fault is 15 km to the southeast), which is poorly constrained on the relatively flat andraised plateau at this location. The GPS dataset is likely to be a good representation of theregional strain. The equal weighted average curve gives a high probability of exceedance of 6
nanostrain/yr (60%), but low probability of exceedance of 8 nanostrain/yr (20%).

Magma intrusion probabilities and volcanic/tectonic strain interactions
The location is approximately 15 km from the nearest active fault and 9 km from the nearest
active volcano (Daibutsu Volcano), on the high plateau linking the Sengan volcanic clusterwith the Japan Sea coast. It is also 15 km from the projected centre of the Sengan volcano
cluster. Short-term, earthquake or eruption scale interactions are less likely to affect thislocation than location 3. Long-term interactions between deep-seated tectonic processes,
such as mantle flow, and future volcanism are likely, but would require an understanding of
which processes are driving uplift of the plateau and its evolution.
The test locations < 30 km from Sengan cluster show the highest hazard, decreasing with
distance from away from the volcano cluster. The inter-quartile range, median, maximum,
logarithmic mean and arithmetic mean are all highest for these localities. The hazards are
only slightly less (55 to 85%) than the hazard at location 3. The alternate dataset analysis forthis particular location (using dataset 3 as a reference dataset), suggests it has fairly highprobabilities for new volcanism in the next 1 Ma, with the logarithmic mean, arithmetic mean
and maximum probabilities being 0.039, 0.026 and 0.210 respectively.

Site investigation needs and siting confidence assessment
This location may be excluded from consideration for a subsequent PIA because it is situatedwithin 8 km of volcanic edifice. However, if further investigations were undertaken, inactive(bedrock) faults have been mapped within 2 km and would need characterisation. Processes
causing plateau uplift need to be understood. The following investigations would be required:detailed site geology investigations (mapping, geomorphology, structural analysis of bedrock);geophysical techniques (GPR, seismics, gravity); determination of rates and mechanisms of
uplift using geochronology, geomorphology, numerical modelling; refinement of tilt and folding
deformation estimates by improved Quaternary mapping and dating of landforms; if potentially
active faults are found then more detailed palaeoseismic techniques (trenching, drilling,surveying) are needed; enhanced GPS network to decrease uncertainties in GPS strain rates;
deployment of microseismic network to obtain detailed distribution of local seismicity including
depths, lineations, rate parameters and focal mechanisms.
Due to the moderate levels of tectonic strain, this location could present moderate
susceptibility to tectonic processes over the next 10,000 years. This location would be
categorised as Level 2: Medium Confidence, requiring detailed assessment of impact
scenarios and a quantitative assessment of risks by the PA team.

- 15 -

Location 5 - GPS
0.4
0.35

Weighted Frequency

0.3
0.25
0.2
0.15
0.1
0.05
0
1

5

9 13 17 21 25 29 33 37 41 45 49 53 57 61 65 69 73 77 81 85 89
nanotrains/ye ar

Location 5 - SEIS
0.4
Weighted Frequency

0.35
0.3
0.25
0.2
0.15
0.1
0.05
1
5
9
13
17
21
25
29
33
37
41
45
49
53
57
61
65
69
73
77
81
85
89

0
nanos trains /ye ar

Location 5 - SURF DEF
0.4

NB. truncatedWeighted Frequency

0.35
0.3
0.25
0.2
0.15
0.1
0.05
0
1 5

9 13 17 21 25 29 33 37 41 45 49 53 57 61 65 69 73 77 81 85 89

nanos trains /yr- 16 -Example Location 6: 39.900N 141.110E
Tectonic setting and geological environment
Example location 6 is 10 km east of the Backbone Range in a wide depression between the
Backbone Range and the Kitakama Block. The location lies 15 km to the east of the volcanic
edifice of Iwate, but 6 km east of Iwate lahars, so would fall within a volcanic exclusion zone.
The nearest active fault is 10 km to the south, and there is also a possibility of faulting along
the west side of the Kitakama Block. The geological formations at the location comprise
Pleistocene volcaniclastic deposits, which are possibly very thick (hundreds of metres?). The
bedrock is uncertain, but may be Miocene andesitic volcanics and possibly Late Palaeozoic to
Mesozoic indurated sediments of the Kitikama Block.
The individual strain indicators are shown at the end of the example description.

Combined strain budgets and probabilities

Location 6

1

Probability of ExceedanceGPS
SURF DEF

0.8

SEIS
EQ W GTD AV

0.6
0.4
0.2
0
1

10

nanostrains/year

100

1000

GPS: High probability of exceedance of 20 nanostrain/yr (60%), but low probability ofexceedance of 50 nanostrain/yr (20%).
SEIS: High probability of exceedance of 4 nanostrain/yr (60%), but low probability ofexceedance of 10 nanostrain/yr (20%).
SURFACE DEFORMATION: High probability of exceedance of 9 nanostrain/yr (60%), but low
probability of exceedance of 12 nanostrain/yr (20%).
STRAIN RATE UNCERTAINTY
Percentile

GPS

Surf Def

Seis

Weight

16th

5

7

3

4

50th

28

8

5

8

84th

50

10

12

29

- 17 -

The equal weighting curve of all three methods provides a measure by which to compare
strain rate models among locations. Seismicity rates may be slightly elevated due to volcanic
processes in that area. The surface deformation dataset is not very discriminating of strain
because of it is dominated by the tilt data (the nearest active fault is 10 km to the south),which is poorly constrained on the relatively flat and raised plateau at this location.As seen in
the GPS histogram there is significant uncertainty because it is closer to the subduction
margin and is probably being influenced by coupling models for the subduction thrust.
However, it is near the eastern margin of the backbone range contractional zone, so the highcontractional strains may be consistent with the tectonic processes. The equal weighted
average curve gives a high probability of exceedance of 8 nanostrain/yr (60%), but low
probability of exceedance of 30 nanostrain/yr (20%).

Magma intrusion probabilities and volcanic/tectonic strain interactions
Location 6 is 20 km east of the main Sengan volcanic cluster centre in a lowland area, on theedge of the high plateau defined by the volcanic front. It is close to Quaternary volcanoes, so
a fairly high hazard is expected and observed. It is only 5 km further from the centre of the
Sengan cluster than location 3, but the hazard is reduced by about a third (Table 2). It also
shows a marked decrease in the range of hazard values produced by the alternative datasets.
Of the 14 locations tested, locality 6 generally has the 4th or 5th highest probabilities relative toother localities. The location is 15 km from the volcanic edifice of Iwate. Stresses induced bymagma intrusion at Iwate volcano, such as eruption or dyking, or deflation events at the
volcano could all impose significant perturbations at this location, such that unknown faults
could be re-activated, or seismicity in the local area might be triggered. The location has a low
but not-insignificant possibility of being inundated by a large lahar event from Iwate volcano.
The example locations <30 km from the Sengan cluster show the highest hazard, decreasing
with distance away from the volcano cluster. The inter-quartile range, median, maximum,
logarithmic mean and arithmetic mean are all highest for these localities. The hazards are
only slightly less (55 to 85%) than the hazard at location 3. Location 6 is ranked fourth out of
the 14 locations in terms of highest probabilities of a new volcano forming in the next 1 Ma.
The alternate dataset analysis for this particular location (using dataset 3 as a reference
dataset), suggests it has fairly high probabilities for new volcanism in the next 1 Ma, with thelogarithmic mean, arithmetic mean and maximum probabilities being 0.037, 0.029 and 0.296
respectively.
There is a long north-south trending fault that is directly in line with this location. Any further
lateral propagation of this fault zone may impact the location.

Site investigation needs and siting confidence assessment
This location may be excluded from consideration for a subsequent PIA because it is situatedwithin the 15 km volcanic exclusion zone. However, if further investigations were undertaken,
the possibility that the nearest active fault may continue northward along the western marginof the Kitakama Block should be investigated. Also, rangefront active faults to the west may
migrate into the basin toward the location. Lahars occur from Iwate volcano 6 km to thesouthwest and should be investigated. The following investigations would be required:detailed site geology investigations (mapping, geomorphology, structural analysis of bedrock);geophysical techniques (GPR, seismics, gravity); determination of rates and mechanisms of
uplift using geochronology, geomorphology, numerical modelling; refinement of tilt and folding
deformation estimates by improved Quaternary mapping and dating of landforms; if potentially
active faults are found then more detailed palaeoseismic techniques (trenching, drilling,surveying) are needed; enhanced GPS network to decrease uncertainties in GPS strain rates;
deployment of microseismic network to obtain detailed distribution of local seismicity including
depths, lineations, rate parameters and focal mechanisms.
Due to the moderate levels of tectonic strain, this location could present moderate
susceptibility to tectonic processes over the next 10,000 years. It would be categorised as
Level 2: Medium Confidence, requiring detailed assessment of impact scenarios and aquantitative assessment of risks by the PA team.

- 18 -

Location 6 - GPS
0.4
0.35

Weighted Frequency

0.3
0.25
0.2
0.15
0.1
0.05
0
1

5

9

13 17 21 25 29 33 37 41 45 49 53 57 61 65 69 73 77 81 85 89
nanostrains/year

Location 6 - SEIS
0.4
Weighted Frequency

0.35
0.3
0.25
0.2
0.15
0.1
0.05
1
5
9
13
17
21
25
29
33
37
41
45
49
53
57
61
65
69
73
77
81
85
89

0
nanos trains /ye ar

Location 6 - SURF DEF
0.4

NB. truncatedWeighted Frequency

0.35
0.3
0.25
0.2
0.15
0.1
0.05
0
1 5

9 13 17 21 25 29 33 37 41 45 49 53 57 61 65 69 73 77 81 85 89

nanos trains /yr- 19 -Example Location 7: 40.000N 139.000E
Tectonic setting and geological environment
Example location 7 is situated offshore, in the Sea of Japan, within the fold and thrust belt
forming an incipient subduction margin. The area is characterised by high historical seismicity
rates and active faults mapped from seismic reflection data. The nearest active fault is 11 km,
but blind thrusts and flexural slip faulting may be very widespread. The nearest Quaternary
active volcano – an isolated submarine cone – is at least 50 km distant. The nearest landbased volcano is Toga Volcano at 60 km to the east. Location 7 is situated within Quaternaryand late Tertiary sedimentary rocks, which could be several km thick.
The seismic strain indicator results are shown below. No GPS or tilt data are available,
although these might be developed from offshore geodetic studies and seismic reflection
profiles, respectively. The surface deformation histogram is truncated (minimum strain <1nanostrain/yr).

Location 7 - SEIS
0.4
Weighted Frequency

0.35
0.3
0.25
0.2
0.15
0.1
0.05
1
5
9
13
17
21
25
29
33
37
41
45
49
53
57
61
65
69
73
77
81
85
89

0
nanos trains /ye ar

Location 7 - SURF DEF
0.4

Weighted Frequency

0.35

NB. truncated

0.3
0.25
0.2
0.15
0.1
0.05
0
1

5

9 13 17 21 25 29 33 37 41 45 49 53 57 61 65 69 73 77 81 85 89

nanos trains /yr- 20 -

Strain budget and probabilities

Location 7
1
Probability of Exceedance

SEIS
SURF DEF

0.8

EQ W GTD AG
0.6
0.4
0.2
0
1

10

100

1000

nanostrains/year
SEIS: High probability of exceedance of 30 nanostrain/yr (60%), but low probability ofexceedance of 50 nanostrain/yr (20%).
SURFACE DEFORMATION: High probability of exceedance of 0.5 nanostrain/yr (60%), but
low probability of exceedance of 1 nanostrain/yr (20%).
STRAIN RATE UNCERTAINTY
Percentile

GPS

Surf Def

Seis

16th

1

1

6

50th

2

1

28

84th

3

1

52

Seismicity rates may be significantly elevated because of recent large earthquakes andaftershock sequences. The surface deformation dataset is limited to active fault data, which ispoorly constrained offshore. The weighted average curve gives a high probability ofexceedance of 1 nanostrain/yr (60%), but low probability of exceedance of 25 nanostrain/yr(20%).

Magma intrusion probabilities and volcanic/tectonic strain interactions
At this location, in the active fault convergence zone on the eastern margin of the Japan Sea,
there is not expected to be any influence of volcanic activity because the nearest Quaternary
active volcano – an isolated submarine cone – is at least 50 km distant. Undersea mappingof volcanic cones may be incomplete, however, and this would need investigation. The
nearest land-based volcano is Toga Volcano, at 60 km to the east.
Localities greater than 70 km away from the Sengan cluster show hazards from 10% to lessthan 1% of the hazard at location 3. These show the lowest maximum and medianprobabilities and also have the smallest ranges. The alternate dataset analysis for thisparticular location (using dataset 3 as a reference dataset), suggests it has very low
probabilities for new volcanism in the next 1 Ma, with the logarithmic mean, arithmetic mean
and maximum probabilities being 0.0016, 0.0005 and 0.014 respectively.

- 21 -

Siting confidence assessment
This location was selected for comparative purposes only. As it is very far offshore,
discussion of site investigation requirements is not appropriate.
Only two datasets are available for evaluation, one of which produces very high strain rates.
Although there are no GPS data available for this offshore location, the high seismic strain
rates are consistent with estimates from GPS that there should be 15-20 mm/yr of
convergence between Northern Honshu and the Japan Sea. Fault map locations areindicative only and require more detailed analysis of seismic reflection data.
Due to the proximity to active faults and high historical seismicity rates, this location could
present high probabilities of some types of tectonic impact over the next 10,000 years, andwould be categorised as Level 3: Limited Confidence.

- 22 -Example Location 8: 40.200N 140.000E
Tectonic setting and geological environment
Example location 8 is located on the Sea of Japan Coast, in the fold and thrust belt
comprising the Sea of Japan contractional zone. The area is marked by high historicalseismicity rates. The nearest active fault is within 3 km and blind thrusts and flexural slip
faulting may be very widespread. The location is north-west of the Sengan Cluster, in the lowlying inter-cluster region, and the nearest volcano (Kampu Volcano) is 32 km distant. Location
8 is situated within Holocene sand, but the underlying bedrock is likely to be Pliocene
sandstone and siltstone.
The individual strain indicators are shown at the end of the example description. The surface
deformation histogram is truncated (maximum strain = 672 nanostrain/yr).

Combined strain budgets and probabilities

Location 8
1

Probability of ExceedanceGPS
SURF DEF

0.8

SEIS
EQ W GTD AV

0.6

0.4

0.2

0
1

10

100

1000

nanos trains /ye ar

GPS: High probability of exceedance of 7 nanostrain/yr (60%), but low probability ofexceedance of 35 nanostrain/yr (20%).
SEIS: High probability of exceedance of 10 nanostrain/yr (60%), but low probability ofexceedance of 35 nanostrain/yr (20%).
SURFACE DEFORMATION: High probability of exceedance of 100 nanostrain/yr (60%), but
low probability of exceedance of 170 nanostrain/yr (20%).
STRAIN RATE UNCERTAINTY
Percentile

GPS

Surf Def

Seis

Weight

16th

2

79

4

4

50th

13

92

14

32

84th

37

158

38

91

- 23 -

The equal weighting curve of all three methods provides a measure by which to compare
strain rate models among locations. As seen in the GPS histogram there is some uncertainty
in the strain rates for this location, due to choice of models for interseismic coupling on
offshore faults in the Sea of Japan convergent zone. This is reflected in the bimodal nature of
the GPS histograms. Seismicity rates are reasonably representative of long-term rates ofseismicity. The surface deformation dataset is dominated by the active fault data, which iswell constrained in this region. Both the fault slip rates and tilt rates are some of the highest in
the Tohoku Region. The equal weighted average curve gives a high probability ofexceedance of 320 nanostrain/yr (60%), but low probability of exceedance of 100
nanostrain/yr (20%).

Magma intrusion probabilities and volcanic/tectonic strain interactions
Volcano-tectonic interactions at this location are unlikely, because the nearest volcano(Kampu Volcano) is 32 km distant and other phenomena that have been linked to volcanism,
such as underplating of magmatic material, leading to an elevated plateau, are lacking.
Location 8 is a coastal site, north-west of the Sengan cluster, in the low-lying inter-cluster
region. There are north-south trending faults nearby.
Localities greater than 70 km away from the Sengan cluster show hazards from 10% to lessthan 1% of the hazard at location 3. These show the lowest maximum and medianprobabilities and also have the smallest ranges. The alternate dataset analysis for thisparticular location (using dataset 3 as a reference dataset), suggests it has fairly low
probabilities for new volcanism in the next 1 Ma, with the logarithmic mean, arithmetic mean
and maximum probabilities being 0.011, 0.0049 and 0.053 respectively.

Site investigation needs and siting confidence assessment
Being within 3 km of the trace of a fault would require careful examination of the accuracy of
the fault mapping, complexity of the faulting and the width of the process zone. An additional
concern is how the fault evolves over the next 10,000 years. If further investigations were
undertaken, the thickness of the poorly consolidated Holocene sediment would need to beconstrained from drilling. This area is an evolving fold and thrust belt associated with anincipient subduction margin which could consist of blind thrusts, tilting, and flexural slip
faulting.
The following investigations would also be required: detailed site geology investigations(mapping, geomorphology, structural analysis of bedrock); geophysical techniques (GPR,
seismics, gravity); determination of rates and mechanisms of uplift using geochronology,
geomorphology, numerical modelling; refinement of tilt and folding deformation estimates byimproved Quaternary mapping and dating of landforms; if potentially active faults are foundthen more detailed palaeoseismic techniques (trenching, drilling, surveying) are needed;enhanced GPS network to decrease uncertainties in GPS strain rates; deployment of
microseismic network to obtain detailed distribution of local seismicity including depths,
lineations, rate parameters and focal mechanisms.
Due to the proximity of an active fault, this location could present high probabilities of some
types of tectonic impact over the next 10,000 years and would be categorised as Level 3:
Limited Confidence.

- 24 -

Location 8 - GPS
0.4
0.35

Weighted Frequency

0.3
0.25
0.2
0.15
0.1
0.05

89

85

81

77

73

69

65

61

57

53

49

45

41

37

33

29

25

21

17

9

13

5

1

0
nanostrains/year

Location 8 - SEIS
0.4

Weighted Frequency

0.35
0.3
0.25
0.2
0.15
0.1
0.05
1
5
9
13
17
21
25
29
33
37
41
45
49
53
57
61
65
69
73
77
81
85
89

0
nanos trains /ye ar

Location 8 - SURF DEF
0.4

Weighted Frequency

0.35

NB. truncated

0.3
0.25
0.2
0.15
0.1
0.05
0
1 5

9 13 17 21 25 29 33 37 41 45 49 53 57 61 65 69 73 77 81 85 89

nanos trains /yr- 25 -Example Location 9: 40.600N 141.500E
Tectonic setting and geological environment
Location 9 is on the Pacific coast, on the east side of the stable forearc (Kitakama) block. It is
one of the closest locations to the subduction thrust. This location is in an area of no
volcanism, and is approximately 95 km from the centre of the Sengan cluster (the densest
volcanic centre in the study area). Hakkoda is the closest volcano, approximately 50 km to thewest-north-west. Location 9 is situated within Holocene gravel, sand, and mud, volcanic ash
and pumice. The bedrock is likely to be Pliocene sandstone and mudstone with conglomerate
and tuff.
The individual strain indicators are shown at the end of the example description. Note that
there are no tilt data, although these might be derived from marine terrace deformation.

Combined strain budgets and probabilities

Location 9

1

Probability of Exceedance

GPS
0.8

SEIS
EQ W GTD AV

0.6
0.4
0.2
0
1

10

nanostrains/year

100

GPS: High probability of exceedance of 5 nanostrain/yr (60%), but low probability ofexceedance of 15 nanostrain/yr (20%).
SEIS: High probability of exceedance of 4 nanostrain/yr (60%), but low probability ofexceedance of 15 nanostrain/yr (20%).
STRAIN RATE UNCERTAINTY
Percentile

GPS

Surf Def

Seis

Weight

16th

3

3

2

2

50th

7

3

7

7

84th

16

3

16

16

- 26 -

1000

The equal weighting curve of all three methods provides a measure by which to compare
strain rate models among locations. As seen in the GPS histogram there is some uncertainty
in the strain rates for this location, due to choice of models for interseismic coupling on the
subduction thrust. Seismicity rates are relatively low, but are probably representative of longterm rates of seismicity. The tilt dataset does not extend to this location, but the presence ofmarine terraces means it could be easily developed and would provide constraints on longterm deformation. The equal weighted average curve gives a high probability of exceedanceof 4.5 nanostrain/yr (60%), but low probability of exceedance of 15 nanostrain/yr (20%).

Magma intrusion probabilities and volcanic/tectonic strain interactions
At this location there is not expected to be any influence of volcanic activity in the arc on
tectonic deformation. The volcanic front has been stable in its location for more than 105 years.
Location 9 is approximately 95 km from the centre of the Sengan cluster (the densest volcanic
centre in the study area). Hakkoda is the closest volcano, approximately 50 km to the westnorth-west.
Localities greater than 70 km away from the Sengan cluster show hazards from 10% to lessthan 1% of the hazard at location 3. These show the lowest maximum and medianprobabilities and also have the smallest ranges. The alternate dataset analysis for thisparticular location (using dataset 3 as a reference dataset), suggests it has fairly low
probabilities for new volcanism in the next 1 Ma, with the logarithmic mean, arithmetic mean
and maximum probabilities being 0.0061, 0.0021 and 0.046 respectively.
This location is fairly close to a mapped off-shore fault which runs parallel to the direction of
the subduction zone. It is likely that there are more off-shore faults which are not mapped
and may interact with this location. It is on the eastern edge of an uplifted plateau where
there is no volcanism and was chosen for investigation to assess the distribution of strain
rates for the so-called ‘stable’ block on the eastern third of Northern Honshu.

Site investigation needs and siting confidence assessment
GPS and seismicity show low rates of strain, well below strain rates of long-term concern. The
mechanism of uplift of marine terraces (possible offshore fault or regional uplift?) needs to beconstrained, and the tilt dataset needs to be extended to this area. The poorly consolidated
Holocene sediments also need to be constrained. This is one of the closest locations to thesubduction thrust, so could be affected by significant ground shaking and tsunami.
The following investigations would also be required: detailed site geology investigations(mapping, geomorphology, structural analysis of bedrock); geophysical techniques (GPR,
seismics, gravity); determination of rates and mechanisms of uplift using geochronology,
geomorphology, numerical modelling; refinement of tilt and folding deformation estimates byimproved Quaternary mapping and dating of landforms; if potentially active faults are foundthen more detailed palaeoseismic techniques (trenching, drilling, surveying) are needed;enhanced GPS network to decrease uncertainties in GPS strain rates; deployment of
microseismic network to obtain detailed distribution of local seismicity including depths,
lineations, rate parameters and focal mechanisms.
Due to the low levels of strain, this location could present low susceptibility to any significant
tectonic hazard over the next 10,000 years, and would be categorised as Level 1: High
Confidence.

- 27 -

Location 9 - GPS
0.4
0.35

Weighted Frequency

0.3
0.25
0.2
0.15
0.1
0.05
0
1

5

9

13 17 21 25 29 33 37 41 45 49 53 57 61 65 69 73 77 81 85 89
nanostrains/year

Location 9 - SEIS
0.4

Weighted frequency

0.35
0.3
0.25
0.2
0.15
0.1
0.05
1
5
9
13
17
21
25
29
33
37
41
45
49
53
57
61
65
69
73
77
81
85
89

0
nanos trains /ye ar

Location 9 - SURF DEF
0.4

NB. truncatedWeighted Frequency

0.35
0.3
0.25
0.2
0.15
0.1
0.05
0
1 5

9 13 17 21 25 29 33 37 41 45 49 53 57 61 65 69 73 77 81 85 89

nanos trains /yr- 28 -Example Location 10: 40.300N 140.500E
Tectonic setting and geological environment
Example location 10 lies between the backbone range and the fold and thrust belt of the Sea
of Japan convergence zone. It is situated on the northern edge of a valley, between twoelevated plateaus. The nearest active fault is 30 km to the ESE and NNE. Several bedrockfaults have been mapped within 0.5 – 5 km. The nearest volcano (Tashiro Volcano) is only 15
km distant. Location 10 is situated on bedrock comprising Miocene lavas and pyroclastic
rocks.
The individual strain indicators are shown at the end of the example description.

Combined strain budgets and probabilities

Location 10
1

Probability of ExceedanceGPS
SURF DEF

0.8

SEIS
EQ W GTD AV

0.6
0.4
0.2
0
1

10

nanostrains/year

100

1000

GPS: High probability of exceedance of 9 nanostrain/yr (60%), but low probability ofexceedance of 15 nanostrain/yr (20%).
SURFACE DEFROMATION: High probability of exceedance of 30 nanostrain/yr (60%), but
low probability of exceedance of 31 nanostrain/yr (20%).
SEIS: High probability of exceedance of 3 nanostrain/yr (60%), but low probability ofexceedance of 6 nanostrain/yr (20%).
STRAIN RATE UNCERTAINTY
Percentile

GPS

Surf Def

Seis

Weight

16th

3

26

2

2

50th

8

29

3

8

84th

16

31

8

29

- 29 -

The equal weighting curve of all three methods provides a measure by which to compare
strain rate models among locations. As seen in the GPS histogram there is some uncertainty
in the strain rates for this location, due to choice of models for interseismic coupling on the
subduction thrust. Seismicity rates are relatively low, but are probably representative of longterm rates of seismicity. The surface deformation dataset is dominated by the tilt data, which
is not well constrained, because it does not quite extend to this location. However, it could be
usefully developed because of the presence of river terraces nearby. The equal weighting
curve (EQ WGTD AV) of the three methods provides a measure by which to compare strain
rate models among sites. The equal weighted average curve gives a high probability ofexceedance of 8 nanostrain/yr (60%), but low probability of exceedance of 30 nanostrain/yr(20%).

Magma intrusion probabilities and volcanic/tectonic strain interactions
Volcano-tectonic interactions at this location are possible because the nearest volcano(Tashiro Volcano) is only 15 km distant. There are numerous bedrock faults close to thelocation. Investigating these features would establish any past volcano-fault interaction, suchas triggered fault slip at the time of major eruptive events at Odate volcano.
Locations that are 30-70 km from the Sengan cluster all show a similar median, range andmaximum, log-mean and arithmetic mean with no tendency for probabilities to change with
distance. The hazard is about 14 to 30% of the hazard at location 3. The alternate dataset
analysis for this particular location (using dataset 3 as a reference dataset), suggests it has
intermediate probabilities for new volcanism in the next 1 Ma, with the logarithmic mean,arithmetic mean and maximum probabilities being 0.021, 0.011 and 0.109 respectively.

Site investigation needs and siting confidence assessment
The GPS and seismicity indicators show low rates of strain, below the strain rates of longterm concern. The tilt dataset needs to be extended to this location and inactive (bedrock)faults need to be investigated.
The following investigations would also be required: detailed site geology investigations(mapping, geomorphology, structural analysis of bedrock); geophysical techniques (GPR,
seismics, gravity); determination of rates and mechanisms of uplift using geochronology,
geomorphology, numerical modelling; refinement of tilt and folding deformation estimates byimproved Quaternary mapping and dating of landforms; if potentially active faults are foundthen more detailed palaeoseismic techniques (trenching, drilling, surveying) are needed;enhanced GPS network to decrease uncertainties in GPS strain rates; deployment of
microseismic network to obtain detailed distribution of local seismicity including depths,
lineations, rate parameters and focal mechanisms.
Due to the uncertainties about the levels of tectonic strain and nearby bedrock faulting, thislocation could present moderate susceptibility to tectonic processes over the next 10,000
years. It would be categorised as Level 2: Medium Confidence, requiring detailed assessmentof impact scenarios and a quantitative assessment of risks by the PA team.

- 30 -

Location 10 - GPS
0.4
0.35

Weighted Frequency

0.3
0.25
0.2
0.15
0.1
0.05
0
1

5

9

13 17 21 25 29 33 37 41 45 49 53 57 61 65 69 73 77 81 85 89
nanostrains/year

Location 10 - SEIS
0.4

Weighted frequency

0.35
0.3
0.25
0.2
0.15
0.1
0.05
1
5
9
13
17
21
25
29
33
37
41
45
49
53
57
61
65
69
73
77
81
85
89

0
nanos trains /ye ar

Location 10 - SURF DEF
0.4

Weighted Frequency

0.35
0.3
0.25
0.2
0.15
0.1
0.05
0
1 5

9 13 17 21 25 29 33 37 41 45 49 53 57 61 65 69 73 77 81 85 89

nanos trains /yr- 31 -Example Location 11: 39.900N 141.000E
Tectonic setting and geological environment
Example location 11 lies on the west side of a wide depression between the backbone range
and the Kitakama Block. The nearest active fault is 9 km to the southwest, but could continue
to the north along the west side of the depression, within 3 km of the location. Location 11 is 6
km north – west of the volcanic edifice of Iwate. It is situated within Pleistocene volcaniclastic
deposits, which are possibly very thick (hundreds of metres?). The bedrock is likely to be
Pleistocene andesite lavas and pyroclastic rocks which are mapped to the north.
The individual strain indicators are shown at the end of the example description. The surface
deformation histogram is truncated (maximum strain = 94 nanostrain/yr).

Combined strain budgets and probabilities

Location 11

1
Probability of ExceedanceGPS
SURF DEF

0.8

SEIS
EQ W GTD AV

0.6
0.4
0.2
0
1

10

100

1000

nanostrains/year
GPS: High probability of exceedance of 9 nanostrain/yr (60%), but low probability ofexceedance of 25 nanostrain/yr (20%).
SEIS: High probability of exceedance of 4 nanostrain/yr (60%), but low probability ofexceedance of 12 nanostrain/yr (20%).
SURFACE DEFORMATION: High probability of exceedance of 10 nanostrain/yr (60%), but
low probability of exceedance of 12 nanostrain/yr (20%).
STRAIN RATE UNCERTAINTY
Percentile

GPS

Surf Def

Seis

Weight

16th

4

8

3

4

50th

10

9

5

9

84th

40

11

14

16

- 32 -

The equal weighting curve of all three methods provides a measure by which to compare
strain rate models among locations. Seismicity rates may be slightly elevated due to volcanic
processes in that area. The surface deformation dataset is dominated by the tilt data, which isprobably representative of long-term deformation at this location. The nearest active fault is 9
km to the south but could extend north to within 3 km. The GPS dataset may be influenced byvolcanic processes, but the rates may be representative of tectonic strain. The equal
weighted average curve gives a high probability of exceedance of 10 nanostrain/yr (60%), but
low probability of exceedance of 12 nanostrain/yr (20%).

Magma intrusion probabilities and volcanic/tectonic strain interactions
This location is 6 km northwest of the volcanic edifice of Iwate. Stresses induced by magma
intrusion at Iwate volcano, by eruption or diking, or deflation events at the volcano could all
impose significant perturbations at the location such that unknown faults could be re-activated,or seismicity in the local area might be triggered. The location lies about 9 km from a mapped
active fault, but could be within 3 km if the fault extends north along the eastern margin of the
Ou-Backbone Range. Near-field (< 10 km) interactions at earthquake or eruptive crisis scale
are likely at this location. Triggered seismicity in the region and on known active faults 4 km
distant are highly likely in the event of major volcanic activity. Conversely, a major earthquake
could trigger a nearby eruption if a magma chamber were susceptible.
The test locations < 30 km from Sengan cluster show the highest hazard, decreasing with
distance from away from the volcano cluster. The inter-quartile range, median, maximum,
logarithmic mean and arithmetic mean are all highest for these localities. The hazards are
only slightly less (55 to 85%) than the hazard at location 3. The alternate dataset analysis forthis particular location (using dataset 3 as a reference dataset), suggests it has fairly highprobabilities for new volcanism in the next 1 Ma, with the logarithmic mean, arithmetic mean
and maximum probabilities being 0.062, 0.040 and 0.683 respectively. This location is thesecond most likely to experience new volcanic interaction according to the alternate dataset
results; 100 times more likely than outlying locations such as 7 and 12.

Site investigation needs and siting confidence assessment
This location may be excluded from consideration for a subsequent PIA because it is justwithin the 15 km volcano exclusion zone. However, if further investigations were undertaken,
the possibility that the nearest active fault may continue northward should be investigated.
Lahars from Iwate volcano 2 km to the south and should also be investigated. The following
investigations would also be required: detailed site geology investigations (mapping,geomorphology, structural analysis of bedrock); geophysical techniques (GPR, seismics,
gravity); determination of rates and mechanisms of uplift using geochronology,
geomorphology, numerical modelling; refinement of tilt and folding deformation estimates byimproved Quaternary mapping and dating of landforms; if potentially active faults are foundthen more detailed palaeoseismic techniques (trenching, drilling, surveying) are needed;enhanced GPS network to decrease uncertainties in GPS strain rates; deployment of
microseismic network to obtain detailed distribution of local seismicity including depths,
lineations, rate parameters and focal mechanisms.
Due to the moderate levels of tectonic strain, this location could present moderate
susceptibility to tectonic processes over the next 10,000 years. This location would be
categorised as Level 2: Medium Confidence, requiring detailed assessment of impact
scenarios and a quantitative assessment of risks by the PA team.

- 33 -

Location 11 - GPS
0.4
0.35

Weighted Frequency

0.3
0.25
0.2
0.15
0.1
0.05

89

85

81

77

73

69

65

61

57

53

49

45

41

37

33

29

25

21

17

9

13

5

1

0
nanostrains/year

Location 11 - SEIS
0.4

Weighted Frequency

0.35
0.3
0.25
0.2
0.15
0.1
0.05
1
5
9
13
17
21
25
29
33
37
41
45
49
53
57
61
65
69
73
77
81
85
89

0
nanos trains /ye ar

Location 11 - SURF DEF
0.4

Weighted Frequency

0.35
0.3
0.25
0.2
0.15

NB. truncated

0.1
0.05
0
1 5

9 13 17 21 25 29 33 37 41 45 49 53 57 61 65 69 73 77 81 85 89

nanos trains /yr- 34 -Example Location 12: 39.650N 142.000E
Tectonic setting and geological environment
Example location 12 is located on the Pacific coast, on the east side of the stable forearc(Kitakama) block. It is the closest of the locations to the subduction thrust. The absence of
Late Pleistocene marine terraces means uplift is not constrained and the embayedmorphology suggests possibility of long-term subsidence. The nearest active fault is 50-80 km
to the west. Location 12 is 100 km from the centre of the Sengan cluster. Iwate Volcano is theclosest Quaternary volcano, located 88 km north-west of Location 12. The location lies within
bedrock comprising Cretaceous and Paleogene intrusives and breccias.
The individual strain indicators are shown at the end of the example description.

Combined strain budgets and probabilities

Location 12

1
Probability of Exceedance

GPS
0.8

SURF DEF
SEIS

0.6

EQ W GTD AV

0.4
0.2
0
1

10

100

1000

nanostrains/year
GPS: High probability of exceedance of 2 nanostrain/yr (60%), but low probability ofexceedance of 7 nanostrain/yr (20%).
SEIS: High probability of exceedance of 8 nanostrain/yr (60%), but low probability ofexceedance of 12 nanostrain/yr (20%).
SURFACE DEFORMATION: High probability of exceedance of 3 nanostrain/yr (60%), but low
probability of exceedance of 14 nanostrain/yr (20%).
STRAIN RATE UNCERTAINTY
Percentile

GPS

Surf Def

Seis

Weight

16th

2

11

4

2

50th

3

12

10

10

84th

7

12

18

13

- 35 -

The equal weighting curve of all three methods provides a measure by which to compare
strain rate models among locations. The GPS dataset is likely to be representative. Seismicity
rates are possibly overestimated due to influence of offshore subduction-related seismicity.
The surface deformation dataset is dominated by the tilt data and is probably representativeof the long-term deformation at this location. The tilting mechanism should be investigated.
The equal weighted average curve gives a high probability of exceedance of 8 nanostrain/yr(60%), but low probability of exceedance of 14 nanostrain/yr (20%).

Magma intrusion probabilities and volcanic/tectonic strain interactions
At this location there is not expected to be any influence of volcanic activity in the arc on
tectonic deformation. The volcanic front has been stable in its location for more than 105 years.
Location 12 is 100 km from the centre of the Sengan cluster. Iwate Volcano is the closest
Quaternary volcano, located 88 km north-west of Location 12.
Localities greater than 70 km away from the Sengan cluster show hazards from 10% to lessthan 1% of the hazard at location 3. These show the lowest maximum and medianprobabilities and also have the smallest ranges. The alternate dataset analysis for thisparticular location (using dataset 3 as a reference dataset), suggests it has very low
probabilities for new volcanism in the next 1 Ma, with the logarithmic mean, arithmetic mean
and maximum probabilities being 0.0020, 0.0004 and 0.024 respectively. This suggests thatlocation 12 is consistently 100 times less likely to develop new volcanic activity in the next
1Ma than locations closer to Sengan. This location was chosen for investigation toinvestigate the distribution of strain rates for the so-called ‘stable’ block on the eastern third ofNorthern Honshu.

Site investigation needs and siting confidence assessment
GPS and seismicity show low rates of strain, well below strain rates of long-term concern. The
long-term tectonics (uplift or subsidence) should be constrained by investigation of marine
terraces or drilling. The mechanism of tilting (regional tectonics versus local faulting) should
also be constrained. This is the closest location to the subduction thrust, so could be affectedby significant ground shaking and tsunami. It is the farthest location from the volcanic front, on
the trench side.
The following investigations would also be required: detailed site geology investigations(mapping, geomorphology, structural analysis of bedrock); geophysical techniques (GPR,
seismics, gravity); determination of rates and mechanisms of uplift using geochronology,
geomorphology, numerical modelling; refinement of tilt and folding deformation estimates byimproved Quaternary mapping and dating of landforms; if potentially active faults are foundthen more detailed palaeoseismic techniques (trenching, drilling, surveying) are needed;enhanced GPS network to decrease uncertainties in GPS strain rates; deployment of
microseismic network to obtain detailed distribution of local seismicity including depths,
lineations, rate parameters and focal mechanisms.
Due to the low levels of strain, this location could present low susceptibility to any significant
tectonic hazard over the next 10,000 years, and would be categorised as Level 1: High
Confidence.

- 36 -

Location 12 - GPS
0.4
0.35

Weighted Frequency

0.3
0.25
0.2
0.15
0.1
0.05

89

85

81

77

73

69

65

61

57

53

49

45

41

37

33

29

25

21

17

9

13

5

1

0
na nostrains/ye ar

Location 12 - SEIS
0.4

Weighted Frequency

0.35
0.3
0.25
0.2
0.15
0.1
0.05
1
5
9
13
17
21
25
29
33
37
41
45
49
53
57
61
65
69
73
77
81
85
89

0
nanos trains /ye ar

Location 12 - SURF DEF
0.4

Weighted Frequency

0.35
0.3
0.25
0.2
0.15
0.1
0.05
0
1 5

9 13 17 21 25 29 33 37 41 45 49 53 57 61 65 69 73 77 81 85 89

nanos trains /yr- 37 -Example Location 13: 38.892N 140.423E
Tectonic setting and geological environment
Example location 13 is situated on the western edge of the backbone range, on the southernedge of an elevated plateau. It is within 3 km of an active fault. Bedrock faults are pervasive
beneath and surrounding the location. Location 13 lies on the southern flank of a plateau
associated with the Chokai volcanic cluster, 16 km from Kobinai Volcano. It is situated on
bedrock comprising Miocene andesite lavas and pyroclastics, acid tuff, sandstone, tuff andconglomerate, granodiorite and quartz diorite.
The individual strain indicators are shown at the end of the example description. The surface
deformation histogram is truncated (maximum strain = 224 nanostrain/yr).

Combined strain budgets and probabilities

Location 13
1
Probability of ExceedanceGPS
SURF DEF

0.8

SEIS
EQ W GTD AV

0.6
0.4
0.2
0
1

10

100

1000

nanostrains/year
GPS: High probability of exceedance of 7 nanostrain/yr (60%), but low probability ofexceedance of 20 nanostrain/yr (20%).
SEIS: High probability of exceedance of 3 nanostrain/yr (60%), but low probability ofexceedance of 9 nanostrain/yr (20%).
SURFACE DEFORMATION: High probability of exceedance of 30 nanostrain/yr (60%), but
low probability of exceedance of 40 nanostrain/yr (20%).
STRAIN RATE UNCERTAINTY
Percentile

GPS

Surf Def

Seis

Weight

th

4

25

3

3

50th

11

35

4

15

th

19

54

9

41

16
84

- 38 -

The equal weighting curve of all three methods provides a measure by which to compare
strain rate models among locations. The GPS dataset may be reflecting strain related tointerseismic coupling processes on nearby faults or may be representing strain occurring at
the location. The dataset is consistent with the location within a contractional strain zone.
Seismicity rates possibly underestimate the long-term rates, give that the location is within 3
km of an active fault, and in the vicinity of other active faults, all of which have long recurrence
intervals. The surface deformation dataset is dominated by active fault data, which is probably
representative. The tilt dataset also shows comparable strain rates at this location. The tilting
mechanism should be investigated. The equal weighted average curve gives a highprobability of exceedance of 8 nanostrain/yr (60%), but low probability of exceedance of 35
nanostrain/yr (20%).

Magma intrusion probabilities and volcanic/tectonic strain interactions
This location lies on the southern flank of a plateau associated with the Chokai volcanic
cluster, 3 km from an active fault, and 16 km from Kobinai Volcano. Near-field (< 10 km)
interactions at earthquake or eruptive crisis scale are possible at this location. Triggered
seismicity in the region and on known active faults 3 km distant are highly likely in the event ofmajor volcanic activity. Conversely, a major earthquake could trigger a nearby eruption if amagma chamber were susceptible. Over time, the interaction of volcanism and tectonism
could promote the migration of either fault activity or the locus of volcanism toward thelocation. Long-term interactions between deep-seated tectonic processes, such as mantle
flow, and future volcanism are likely, but would require an understanding of which processes
are driving uplift of the plateau and its evolution.
Localities greater than 70 km away from the Sengan cluster show hazards from 10% to lessthan 1% of the hazard at location 3. These show the lowest maximum and medianprobabilities and also have the smallest ranges. However, location 13 is an exception withelevated hazard (30% compared to location 3). The alternate dataset analysis for thisparticular location (using dataset 3 as a reference dataset), suggests it has fairly highintermediate probabilities for new volcanism in the next 1 Ma, with the logarithmic mean,arithmetic mean and maximum probabilities being 0.024, 0.014 and 0.154 respectively.

Site investigation needs and siting confidence assessment
Being within 3 km of the trace of the fault would require careful examination of the accuracy of
the fault mapping, complexity of the faulting and the width of the process zone. An additional
concern is how the fault evolves over the next 10,000 years. If further investigations were
undertaken, bedrock faults should be examined. The sandstone and mudstone bedrock may
also provide evidence of long-term brittle deformation. A magnitude-frequency analysis ofseismicity would test whether historical seismicity adequately represents the long term large
earthquake rates.
The following investigations would also be required: detailed site geology investigations(mapping, geomorphology, structural analysis of bedrock); geophysical techniques (GPR,
seismics, gravity); determination of rates and mechanisms of uplift using geochronology,
geomorphology, numerical modelling; refinement of tilt and folding deformation estimates byimproved Quaternary mapping and dating of landforms; if potentially active faults are foundthen more detailed palaeoseismic techniques (trenching, drilling, surveying) are needed;enhanced GPS network to decrease uncertainties in GPS strain rates; deployment of
microseismic network to obtain detailed distribution of local seismicity including depths,
lineations, rate parameters and focal mechanisms.
Due to the proximity to an active fault, this location could present high probabilities of some
types of tectonic impact over the next 10,000 years, and would be categorised as Level 3:
Limited Confidence.

- 39 -

Location 13 - GPS
0.4
0.35

Weighted Frequency

0.3
0.25
0.2
0.15
0.1
0.05
0
1

5

9

13 17 21 25 29 33 37 41 45 49 53 57 61 65 69 73 77 81 85 89
nanostrains/year

Location 13 - SEIS
0.4

Weighted Frequency

0.35
0.3
0.25
0.2
0.15
0.1
0.05
1
5
9
13
17
21
25
29
33
37
41
45
49
53
57
61
65
69
73
77
81
85
89

0
nanos trains /ye ar

Location 13 - SURF DEF
0.4

Weighted Frequency

0.35
0.3
0.25
0.2
0.15
0.1
0.05
0
1

5

9 13 17 21 25 29 33 37 41 45 49 53 57 61 65 69 73 77 81 85 89

nanostrains/yr- 40 -Example Location 14: 39.387N 140.303E
Tectonic setting and geological environment
Example location 14 lies on the eastern flank of an uplifted range in the fold and thrust belt
related to the Sea of Japan convergence zone. The location is in a gap region between
volcanic clusters. The closest active faults are 5-6 km to the east and south, but they may
extend further north, close to the location. Bedrock faults have been mapped within 1.5 km.
The nearest volcano (Chokai Volcano) is 35 km to the south-west. Location 14 is situatedwithin bedrock that comprises Late Miocene shale and tuff, mudstone and sandstone.
The individual strain indicators are shown at the end of the example description.

Combined strain budgets and probabilities

Location 14
1

Probability of Exceedance

GPS
0.8

SURF DEF
SEIS

0.6

EQ W GTD AV

0.4
0.2
0
1

10

nanostrains/year

100

1000

GPS: High probability of exceedance of 4 nanostrain/yr (60%), but low probability ofexceedance of 12 nanostrain/yr (20%).
SEIS: High probability of exceedance of 3 nanostrain/yr (60%), but low probability ofexceedance of 8 nanostrain/yr (20%).
SURFACE DEFORMATION: High probability of exceedance of 20 nanostrain/yr (60%), but
low probability of exceedance of 22 nanostrain/yr (20%).
STRAIN RATE UNCERTAINTY
Percentile

GPS

Surf Def

Seis

Weight

16th

2

18

3

3

50th

6

20

4

8

84th

13

21

9

20

- 41 -

The equal weighting curve of all three methods provides a measure by which to compare
strain rate models among locations. The GPS dataset may be reflecting strain related tointerseismic coupling processes on nearby faults or may be representing strain occurring at
the location. The dataset is consistent with the location within a contractional strain zone.
Seismicity rates possibly underestimate the long term rates, give that the location is within 5
km of an active fault, and in the vicinity of other active faults, all of which have long recurrence
intervals. The surface deformation dataset is dominated by the tilt data, which is not wellconstrained due to absence of river terraces. However, it is consistent with presence of active
faults, so is probably representative. The equal weighted average curve gives a highprobability of exceedance of 5 nanostrain/yr (60%), but low probability of exceedance of 20
nanostrain/yr (20%).

Magma intrusion probabilities and volcanic/tectonic strain interactions
Volcano-tectonic interactions at this location are unlikely because the nearest volcano(Chokai Volcano) is 35 km to the south-west and other phenomena that have been linked tovolcanism, such underplating of magmatic material leading to an elevated plateau, are lacking.
Localities greater than 70 km away from the Sengan cluster show hazards from 10% to lessthan 1% of the hazard at location 3. These show the lowest maximum and medianprobabilities and also have the smallest ranges. The alternate dataset analysis for thisparticular location (using dataset 3 as a reference dataset), suggests it has fairly low
probabilities for new volcanism in the next 1 Ma, with the logarithmic mean, arithmetic mean
and maximum probabilities being 0.013, 0.0048 and 0.081 respectively.

Site investigation needs and siting confidence assessment
The possible northward extension of the active faults and the bedrock faults should be
investigated. The presence of Tertiary sediments may provide evidence of long-term brittle
deformation. A magnitude-frequency analysis of seismicity would test whether historicalseismicity adequately represents the long term large earthquake rates.
The following investigations would also be required: detailed site geology investigations(mapping, geomorphology, structural analysis of bedrock); geophysical techniques (GPR,
seismics, gravity); determination of rates and mechanisms of uplift using geochronology,
geomorphology, numerical modelling; refinement of tilt and folding deformation estimates byimproved Quaternary mapping and dating of landforms; if potentially active faults are foundthen more detailed palaeoseismic techniques (trenching, drilling, surveying) are needed;enhanced GPS network to decrease uncertainties in GPS strain rates; deployment of
microseismic network to obtain detailed distribution of local seismicity including depths,
lineations, rate parameters and focal mechanisms.
Due to the moderate levels of tectonic strain, this location could present moderate
susceptibility to tectonic processes over the next 10,000 years. This location would be
categorised as Level 2: Medium Confidence, requiring detailed assessment of impact
scenarios and a quantitative assessment of risks by the PA team.

- 42 -

Location 14 - GPS
0.4
0.35

Weighted Frequency

0.3
0.25
0.2
0.15
0.1
0.05
0
1

5

9 13 17 21 25 29 33 37 41 45 49 53 57 61 65 69 73 77 81 85 89
nanostrains/year

Location 14 - SEIS
0.4

Weighted Frequency

0.35
0.3
0.25
0.2
0.15
0.1
0.05
1
5
9
13
17
21
25
29
33
37
41
45
49
53
57
61
65
69
73
77
81
85
89

0
nanos trains /ye ar

Location 14 - SURF DEF
0.4

Weighted Frequency

0.35
0.3
0.25
0.2
0.15
0.1
0.05
0
1 5

9 13 17 21 25 29 33 37 41 45 49 53 57 61 65 69 73 77 81 85 89

nanos trains /yr- 43 -