CRUSTAL STRUCTURE OF THE MANGYSTAU REGION, WESTERN KAZAKHSTAN Luis Bernardo Martinetti By A THESIS Submitted to Michigan State University in partial fulfillment of the requirements for the degree of Geological Sciences – Master of Science 2019 ABSTRACT By Luis Bernardo Martinetti CRUSTAL STRUCTURE OF THE MANGYSTAU REGION, WESTERN KAZAKHSTAN The Mangystau Region, in southwestern Kazakhstan, is a region of complex structure resulting from terrane collision and rifting in the Late Paleozoic to Early Triassic. The crustal and upper mantle structure was intensely studied with seismic methods during the 1960’s, producing competing interpretations. In 2016, Michigan State University and the Institute of Geophysical Research of Kazakhstan installed 10 broadband seismic stations throughout the region for about 20 days. With the data collected, receiver functions were produced and inverted jointly with surface wave dispersion curves to create shear wave velocity models that estimate the crustal structure. The results are interpreted as a sedimentary basin that is 3 km in North Ustyurt and Central Mangyshlak, and thickening to 5 km in South Mangyshlak. The Permian-Triassic Complex is found throughout all the Mangystau region varying in thickness; 7.5 km in North Ustyurt and Central Mangyshlak, and 4 km in South Mangyshlak. The granitic upper basement thickens from 10 km in North Ustyurt to 13 km in Central Mangyshlak and South Mangyshlak. Since the Moho is not clear in this study, the basaltic lower basement thickness is unknown. Results are in general agreement with prior studies, although some differences exist, perhaps induced by the smoothing parameter during the joint inversion causing sharp boundaries to be smeared in the results. More seismic data is necessary to study this region thoroughly with this inversion method. ACKNOWLEDGEMENTS I would like to thank Dr. Kevin Mackey for all his guidance, constructive feedback, support, and friendship. Second, I would like to thank my family, Luis Kike, and my friends for supporting me when I would encounter hard times. A very special thanks to Dr. Kazuya Fujita for teaching me how to interpret Russian literature and having fruitful discussion of the geology and tectonics. Thanks to Dr. Allen McNamara for serving as a member of my committee. I would like to thank Cristo Ramirez and the group from Penn State for teaching me the receiver function methodology. A great thanks to Chengping Chai and Monica Maceira for inviting me to the University of Tennessee and teaching me how to apply the joint inversions as well as very patiently walking me though the process. I would also like to thank Vaclav Kuna for helping and teaching me how to produce the maps of the region using Generic Mapping Tools [GMT] (Wessel et al., 2013). This deployment of stations was possible with the cooperation of the Institute of Geophysical Research of Kazakhstan, with station deployment funding from the U.S. Department of State, Bureau of Arms Control, Verification, and Compliance. The financial assistantships funding my master’s degree were provided by The Graduate School of Michigan State University and the Michigan State University Department of Earth and Environmental Sciences, and the Michigan Space Grant Consortium. I really want to thank Dr. Michael Velbel for taking me as an undergraduate research assistant during the summer of 2014 and exposing me to Michigan State; thank you for believing in me and guiding me to what my real passion is. Above everything, I would like to thank my mom for being my emotional support and for always believing in me. iii TABLE OF CONTENTS LIST OF TABLES.........................................................................................................................................................v LIST OF FIGURES......................................................................................................................................................vi INTRODUCTION.........................................................................................................................................................1 GEOLOGIC SETTING.................................................................................................................................................3 North Ustyurt.........................................................................................................................................................3 Central Mangyshlak............................................................................................................................................4 South Mangyshlak...............................................................................................................................................5 Tectonic Models.....................................................................................................................................................5 PREVIOUS CRUSTAL STRUCTURE STUDIES.................................................................................................9 Crustal structure of North Ustyurt............................................................................................................11 Crustal structure of Central Mangyshlak................................................................................................15 Crustal structure of South Mangyshlak...................................................................................................17 DATA COLLECTION AND METHODOLOGY .................................................................................................19 Seismic station deployment and data collection.................................................................................19 Methodologies....................................................................................................................................................21 Receiver functions.......................................................................................................................................21 Processing..................................................................................................................................................23 Time Domain Iterative Deconvolution..........................................................................................25 Stacking.......................................................................................................................................................26 Surface Wave dispersion..........................................................................................................................28 Initial shear wave velocity structure models..................................................................................33 Joint inversions of receiver functions and surface wave dispersion.....................................33 RESULTS FOR RECEIVER FUNCTIONS AND LINEARIZED JOINT INVERSION OF RECEIVER FUNCTIONS AND SURFACE WAVE DISPERSION......................................................................................40 Receiver functions results.............................................................................................................................40 Joint inversion of receiver functions with surface wave dispersion results...........................59 Qualitative effect of data quality on results...........................................................................................68 DISCUSSION OF RESULTS...................................................................................................................................72 North Ustyurt......................................................................................................................................................72 Central Mangyshlak..........................................................................................................................................77 South Mangyshlak.............................................................................................................................................79 Moho discontinuity..........................................................................................................................................81 CONCLUSIONS..........................................................................................................................................................83 BIBLIOGRAPHY.......................................................................................................................................................86 iv LIST OF TABLES Table 1: Seismic stations used in this study...............................................................................................20 Table 2: Earthquakes used for receiver function estimations............................................................22 v LIST OF FIGURES Figure 1: Shaded relief map of the Mangystau Region, southwest Kazakhstan, with seismic stations..........................................................................................................................................................................2 Figure 2: Plate tectonic evolution cross-sections of the region from the Carboniferous to Cretaceous....................................................................................................................................................................6 Figure 3: Strike slip shearing and stacking of the Silk Road Arc..........................................................8 Figure 4: Cross sections from seismic refraction studies in the Mangystau Region.................10 Figure 5: Map of the topography of the Mohorovičić discontinuity.................................................12 Figure 6: Interpretation of wide-angle seismic reflection profile of the Aral Sea......................14 Figure 7: Cross-section showing relative ages of the geologic units in the region....................16 Figure 8: Cross section showing folding in the western section of Central Mangyshlak from seismic imaging.......................................................................................................................................................18 Figure 9: Idealized ray diagram and receiver function of a teleseismic earthquake................24 Figure 10: Example of receiver functions with different Gaussian filters from earthquake detected in station WKZ10.................................................................................................................................27 Figure 11: Shear wave sensitivity kernels for Love and Rayleigh waves......................................30 Figure 12: Surface wave dispersion curves for all the seismic stations.........................................32 Figure 13: Diagram showing the range of “a priori” shear wave velocity structure models used..............................................................................................................................................................................34 Figure 14: Comparison of different smoothness factors from 0.15 to 5.00..................................38 Figure 15: Individual receiver functions from all the stations...........................................................41 Figure 16: Individual and stacked receiver functions of station WKZ09.......................................43 Figure 17: Individual and stacked receiver functions of station WKZ10.......................................44 Figure 18: Stack from group of stations WKZ09 and WKZ10.............................................................45 Figure 19: Individual and stacked receiver functions of station WKZ01.......................................47 vi Figure 20: Individual and stacked receiver functions of station WKZ02.......................................48 Figure 21: Stack from group of stations WKZ01 and WKZ02.............................................................49 Figure 22: Individual and stacked receiver functions of station WKZ03.......................................50 Figure 23: Receiver function of station WKZ04........................................................................................52 Figure 24: Stack from group of stations WKZ03 and WKZ04.............................................................53 Figure 25: Individual and stacked receiver functions of station WKZ05.......................................54 Figure 26: Individual and stacked receiver functions of station WKZ06.......................................55 Figure 27: Stack from group of stations WKZ05 and WKZ06.............................................................57 Figure 28: Individual and stacked receiver functions of station WKZ07.......................................58 Figure 29: Joint inversion result for group of stations WKZ09 and WKZ10.................................62 Figure 30: Joint inversion result for group of stations WKZ01 and WKZ02.................................63 Figure 31: Joint inversion result for group of stations WKZ03 and WKZ04.................................65 Figure 32: Joint inversion result for group of stations WKZ05 and WKZ06.................................67 Figure 33: Joint inversion result for station WKZ07...............................................................................69 Figure 34: Columns representing the crustal structure of the Mangystau Region from results..........................................................................................................................................................................73 Figure 35: Comparison of a crustal column from our results in North Ustyurt and a cross- section from Sheikh-Zade (1996)....................................................................................................................74 Figure 36: Comparison of crustal column from our results of North Ustyurt and Dimakov (1968) Profile I........................................................................................................................................................75 Figure 37: Comparison of crustal column from our results of Central Mangyshlak and Dimakov (1968) Profile I....................................................................................................................................78 Figure 38: Comparison of crustal columns from our results of South Mangyshlak and Dimakov (1968) Profile I....................................................................................................................................80 vii INTRODUCTION The purpose of this study is to investigate the variations of the crustal structure across the Mangystau Region (Figure 1) in southwestern Kazakhstan by computing receiver functions from teleseismic earthquakes and inverting them jointly with surface wave dispersion. From August through September of 2016, Michigan State University and the Institute of Geophysical Research of Kazakhstan carried out a temporary deployment of 10 broadband seismic stations for 20 days with the goal of quantifying seismic noise and event detection thresholds throughout the region. The area of the network deployment is approximately bounded by 41-45°N and 51-56°E. This area encompasses part of the Mangystau Province, between the Caspian and Aral Seas. The improved determination of the crustal structure of this region can contribute to numerous objectives. These include understanding the regional tectonic setting, for which there is some disagreement, determine the S-wave velocity structure of the region which has never been directly investigated (I. Sokolova, Pers. Com.), and improve monitoring of seismic events in the region. The knowledge of the crustal structure gained from this region can also help identify any correlations between the structures of this region and the Scythian Platform across the Caspian Sea, where a major effort to understand the young Caucasus Mountains to its south is ongoing. 1 Figure 1: Shaded relief map of the Mangystau Region, southwest Kazakhstan, with seismic stations. The red triangles are the locations of the seismic stations deployed for this study. (top) The North Mangyshlak Fault is the boundary between the Turan Platform (south) and Eurasia Plate (north). The fault also acts as the northern boundary between Central Mangyshlak and North Ustyurt. The Zhetybay-Uzen’ Fault is the boundary between Central Mangyshlak and South Mangyshlak (bottom). The inset reference map showing a broader region. The blue area is the Mangystau Region, and the red box is the location of the top map. 2 GEOLOGIC SETTING The Mangystau Region is located on the southern edge of the present day Eurasian Plate, and includes a broad deformation zone that forms part of the Alpine-Himalayan fold- belt system (Bird, 2003). This region has normally been considered part of the Turan Platform, which extends to the southeastern edge of the East European Platform to the north, the Kopet-Dag mountains in Turkmenistan to the south, bounded by the Caspian Sea to the West, and extending to the Pamir mountains to the east (Khain, 1994). However, different authors have used different naming conventions for geologic and tectonic structures and locations, thus care is needed. This region can be divided into three major, fault bounded, sub- regions (Popkov, 1986): the North Ustyurt region to the north, the South Mangyshlak region to the south, and the Central Mangyshlak (also called Central Ustyurt) region in between the other two. The North Ustyurt and Central Mangyshlak are separated by the North Mangyshlak Fault (also called Ustyurt or North Karatau fault), although its location differs some between authors (Dimakov, 1968; Mstislavskiy et al., 1980; Thomas et al., 1999). Central Mangyshlak and South Mangyshlak are separated by the Zhetybay-Uzen’ fault, which is parallel to the North Mangyshlak Fault (Yuferov et al., 1973 and Mstislavskiy et al., 1980). The basic geologic characteristics from prior studies are summarized below. The basement in this region is presumed to be Pre-Cambrian in age, based on data from dated rocks to the north (Mstislavskiy et al., 1980; Khain, 1994). This basement is overlain by a deformed and metamorphosed Middle Paleozoic basement, which might be part of the East European Platform (Khain ,1994). Boreholes studies by Popkov (1995) 3 North Ustyurt Central Mangyshlak show Carboniferous carbonate-terrigenous-volcanic deposits in some places and only carbonate on others, overlain by Late Permian to Early Triassic redbeds interbedded with tuffs and lavas of mafic to intermediate composition. The shallower and younger part of North Ustyurt is formed by the succession of Lower to Middle Jurassic terrigenous rocks, Late Jurassic marine carbonate and clastic rocks, Early Cretaceous redbeds, mid-Cretaceous marine clastics, Late Cretaceous to Early Paleogene marine carbonates and clastics, and Oligocene clayey rocks (Khain, 1994). The basement of the Central Mangyshlak region does not differ from the North Ustyurt. However, Central Mangyshlak can be divided into a north and south sections separated by the South Karatau Fault (Dimakov, 1968; Mstislavskiy et al., 1980). The northern section of Central Mangyshlak is uplifted and exposes highly deformed and overturned Permian and Triassic rocks. These rocks are thought to have formed during high rates of sedimentation forming marine to alluvial volcanoclastic fans. Some authors think that these deposits might have been formed in a back-arc basin regime (Thomas et al.,1999). The rocks were deformed, faulted, and uplifted in the Late Triassic due to compressional deformation (Zonenshain et al., 1990; Khain, 1994; Thomas et al., 1999). From the Jurassic to Quaternary, the uplifted and Central Mangyshlak has been providing sedimentation for other surrounding areas. The Permian and Triassic rocks are present in the south section of Central Mangyshlak and extend to the Zhetybay-Uzen’ Fault, but they are not as deformed. The Late Permian to Early Triassic basin here is overlain by Cretaceous to Quaternary sedimentary rocks (Mstislavskiy et al., 1980). 4 South Mangyshlak The northern part of South Mangyshlak is inferred to have the Precambrian basement under a series of volcanic arches and basins. The upper crust consists of Late Permian to Early Triassic coarse clastic rocks, thinning southward to essentially non- existent (Thomas et al., 1999). Overlaying it, are terrigenous Early to Middle Jurassic deposit followed by lagoonal-shallow marine Late Jurassic deposits, that also pinch out to the south. It is important to note that the mid-Cretaceous to Oligocene formations (terrigenous and marine deposit alternations) found in North Ustyurt are also found in this region (Khain, 1994). Jurassic and younger deposits seem to have varying degrees of deformation (Khain, 1994). Also, from borehole data, the Permian-Triassic rocks found here are similar in stratigraphy to those found in the northern part of Central Mangyshlak (Natal’in and Şengör, 2005). The southern part of South Mangyshlak is also referred to as the “Karabogaz region”. This is where the sedimentary basins thin out, and base of the crust is found higher. This region includes pre-Permian to Cretaceous and younger formations, similar to elsewhere. Two general models have been proposed to explain the tectonic evolution of the region. The model from Zonenshain et al. (1990), Gaetani et al. (1998), and Golonka (2004) suggest that the Permian-Triassic Complex rocks in Central Mangyshlak were deposited in an extensional back-arc basin, caused by subduction of the Paleotethys under the Turan Platform in the Late Triassic, and forming the Kopet Dag mountains at the same time (Figure 2). This rift separated North Ustyurt from South Mangyshlak, depositing most of the Permian-Triassic sediments in Central Mangyshlak. Later in the Late Triassic, a series of Tectonic models 5 Figure 2: Plate tectonic evolution cross-sections of the region from the Carboniferous to Cretaceous. Each cross section represents the movement of plates at a specific time period and epoch. The left side of the cross-sections is present-day south, and the right is North. Adapted from Golonka (2004). 6 microplates collided to the south of the Turan Platform, closing the Paleotethys along the Kopet Dag Mountains. This compression caused the Central Mangyshlak basin to close. Later accretion of microcontinents on the south edge of the Platform added to the compressional stress of the region causing the deformation and uplifting of the Central Mangyslak in the Late Jurassic. The other proposed model, from Natal’in and Şengör (2005), is based on the idea that the series of alternating volcanic arcs and basins on the margin of the Turan Platform (Central Mangyshlak and South Mangyshlak) are fragments from the Silk Road Arc. This was a magmatic volcanic arc that extended from northern China to the Black Sea, caused by the subduction of the Paleotethys Ocean. Due to oblique convergence of the ocean, the arc was fragmented by strike-slip faults, translated along the southern edge of Eurasia, and accreted to the present North Ustyurt forming this alternating arc-basin pattern in the Late Triassic (Figure 3). Both models agree in that the northern Central Mangyshlak was uplifted due to compressional forces in the Late Triassic, and that the basement in North Ustyurt is similar in many aspects to that in Central Mangyshlak and South Mangyshlak. In terms of crustal structure, each model would suggest something different about South Mangyshlak: the back-arc model suggests that South Mangyshlak and North Ustyurt have same similar crustal structure since these regions separated because of the rifting, and the Silk-Road Arc restacking suggests that the crustal structure in South Mangyshlak should be different than the North Ustyurt, due to the alternating arcs and basins that came from elsewhere. 7 Figure 3: Strike slip shearing and stacking of the Silk Road Arc. A map-view evolution of stacking of arcs and forearcs as suggested by Natal’in and Şengör (2005). The forearc and arc blocks are from a single arc instead of multiple microplates colliding. The numbers to the right indicate sequential periods of development of the evolution of the faulting creating the alternating pattern. From Natal’in and Şengör (2005). 8 PREVIOUS CRUSTAL STRUCTURE STUDIES During the 1960’s, Soviet scientists studied the crust and upper mantle structure using Deep Seismic Sounding (DSS) throughout the Soviet Union (e.g., Ryaboy, 1989). This method focuses on mapping the overall crustal structure using refracted waves and images the Moho and upper mantle. These DSS profiles extending for thousands of kilometers were collected using earthquakes and Peaceful Nuclear Explosions (PNEs) as sources. The first profile was conducted in this study area in 1964-63 and extended for 625 km from the Kopet-Dag to the Aral Sea (Ryaboy 1967). Additional profiles were collected in the Mangystau Region directly crossing our study region that determined seismic velocities, crustal structure, and Moho depth (Volvovsky et al., 1966; Dimakov, 1968). The Dimakov (1968) study provides the majority of the geophysical information and geologic interpretations made in the Mangystau Region from two seismic profiles (Figure 4). The two profiles include seismic velocities, layer thickness, and structure for the western part of the Mangystau region. These profiles image the crust where our stations WKZ03, WKZ05, and WKZ04 were deployed. Dimakov (1968) provides P-wave velocities in his seismic profiles, but they represent mean velocity from the surface to the bottom of the layer, averaging the upper velocities. He also provides the individual layer velocities, which are the ones we use in this study: 3.05 km/s for upper sedimentary layer, 6.00 km/s for granitic-metamorphic layer (our upper basement), 6.80 km/s for the basaltic layer (our lower basement), and 8.20km/s for upper mantle. However, the velocity for the Permian- Triassic Complex was not provided, thus we calculated the layer velocity from the mean velocity (surface to base of layer) to have velocities of 4.98 – 5.20 km/s. Mstislavskiy et al. (1980) mention that the layer velocity of the Permian-Triassic Complex is 5.00 to 5.20 9 Figure 4: Cross sections from seismic refraction studies in the Mangystau Region. Two cross sections modified from the Dimakov (1968) study where two seismic lines crossed and imaged through South Mangyshlak and Central Mangyshlak. The cross sections are shown on the map as purple lines. Profile I is on top and Profile II is on bottom. The yellow boxes on Profile I are the sections used to compare with the columns. The Vp represents the average P-wave velocity for each layer, contrary from the original figure from Dimakov (1968) which shows the average velocity from the surface to the base of the layer. 10 km/s. Since both studies measured velocities in the same range, we used a layer velocity of 5.10 km/s. For this dataset, we use the standard Vp/Vs of 1.70 to keep consistency with the rest of this study. The map of crustal thickness (Moho depth) produced by Volvovsky et al. (1966) shows variation across the region using DSS profiles that extends from near the border of Turkmenistan and Iran to the northeast extent of the Caspian Sea (Figure 5). A profile that crosses our study area in Figure 5 shows a Moho about 40 km deep throughout all three sub-regions, while it sharply decreases to the south. Other studies (Babadjanov et al., 1989, 1991; Sheikh-Zade, 1996; Abetov et al., 2015) suggest that there is a gradational Moho boundary instead of a distinct one throughout all of the Mangystau region. In the next sections, the previous studies on the crustal structure of each sub-region (North Ustyurt, Central Mangyshlak, and South Mangyshlak) are summarized. For the purpose of this summary, we will focus on the crustal structure, the basement, variations within the Permian-Triassic sedimentary Complex, and the Moho. Information on the composition, seismic velocity, and layer thickness are available for most of the areas. One has to be careful interpreting the mean Vp velocities from Russian literature, as they represent the average velocity from the surface to the base of the layer, incorporating layers with different velocities. All the velocities shown here are layer velocities, meaning that they only represent the seismic velocity of the mentioned layer. Two seismic studies were conducted in North Ustyurt, but neither of these are near where our stations were deployed. The far western part of the basin is imaged by both seismic profiles (Figure 4) form Dimakov (1968), while the eastern part is covered by a 200 Crustal structure of North Ustyurt 11 Figure 5: Map of the topography of the Mohorovičić discontinuity. The red numbers represent the Moho depths and the contour intervals are 10 km, and the red box bounds our study area. The DSS line of interest is the one crossing the red box, with absolute depths along the profile. Adapted from Volvovsky et al. (1966). 12 km long SW-NE seismic profile (named Aral-1) conducted by Babadjanov et al. (1989) in the Aral Sea in 1978, and later reinterpreted by Sheikh-Zade (1996) (Figure 6). Unfortunately, the Aral-1 profile is distant away from our stations and may not be representative of our study area. Based on Profile II from Dimakov (1968), since it is the closest to our stations, the top of the North Ustyurt begins with a 2 km deep Jurassic to Quaternary sedimentary basin with layer Vp=3.05 km/s (Vs=1.80 km/s). The underlying unit is a 5 km thick Permian- Triassic Complex with layer Vp=5.10 km/s (Vs=3.00 km/s). A shallow fault is present and both the younger basin and Permian-Triassic Complex increase 1 km in thickness. The upper basement begins at 7 km depth extending to about 27 km, with a total thickness of 20 km and a layer Vp=6.00 km/s (Vs=3.53 km/s). The lower basement extends from 27 km to the base of the crust at 43 km, with a total thickness of 16 km and a layer Vp=6.80 km/s (Vs=4.00 km/s). The Moho is located at a depth of 43 km and it has a velocity of Vp=8.20 km/s. The velocities from the Aral-1 profile (Figure 6) are all also P-wave layer velocities at depth. This profile shows an overall 10 km thick sedimentary cover with two layers: a 3 km upper layer with Vp=3.15 km/s and lower 7 km layer with Vp=5.90-6.50 km/s (Vs=1.85 km/s and 3.47-3.82 km/s) velocities, overlaying a basement with Vp velocities of 6.40-6.60 km/s. The basement layer velocities range from (top to bottom) Vp=6.40-6.6 km/s (Vs=3.76-3.88 km/s) in the upper ~15 km, to Vp=6.30-7.10 km/s (Vs=3.71-4.18 km/s) in the middle 7 km, to about a Vp=6.7-7.5 km/s (Vs=3.94-4.41 km/s) in the lower 5 km of the crust in the Moho transition zone (Sheikh-Zade, 1996), with an overall crustal thickness of 13 Figure 6: Interpretation of wide-angle seismic reflection profile of the Aral Sea. The seismic profile was conducted in the Aral Sea roughly from southwest to northeast (yellow line in the map). The NE part of the profile shows a deep basin (~10 km). Mt and Mb represent the top and bottom boundaries of the Moho transition zone (~ 5 km thick), and the numbers represent layer Vp in km/s at the plotted depth. Cross-section adapted from Sheikh-Zade (1996). Xcv 14 Crustal structure of Central Mangyshlak 35-37 km. Note that there are lateral variation of velocities that mainly increase from SW to NE as the geology changes northward away from Central Mangyshlak. There are some differences between the information provided by the Aral-1 profile and the Dimakov (1968) Profile II. The layer velocities, layer thicknesses, and structure are very similar throughout the whole crust for both models, except for the crustal thickness. A sharp Moho is mapped at 43 km in depth by Dimakov (1968) while a gradational Moho was mapped at 35 – 37 km by Sheikh-Zade (1996). The upper crust in Central Mangyshlak and South Mangyshlak, are dominated by sedimentary layers ranging from the Palaeozoic to the Quaternary ages and have an overall thickness up to 10 km (Thomas et al., 1999), cut by a series of faults (Popkov, 1986) (Figure 7). Due to the complex faulting and folding of the region, the sedimentary cover is not flat throughout, consistent with gravity studies from Volvovsky et al. (1966), shown by Mstislavskiy et al. (1980). The crustal structure previously determined in Central Mangyshlak from Dimakov (1968) near the location of where our station WKZ03 was installed, begins from the surface with 3 km Jurassic to Quaternary sedimentary layer with a layer Vp=3.05 km/s (Vs=1.79 km/s). The folded Permian-Triassic Complex has a thickness of about 12 km with a layer Vp=5.10 km/s (Vs=3.00 km/s) and is observed as slightly metamorphosed sedimentary rocks at the surface. The upper basement section has a thickness of 9 to 10 km with a layer Vp=6.00 km/s (Vs=3.53 km/s). The lower crustal section has a thickness of about 16 km with a layer Vp=6.80 km/s (Vs=4.00 km/s). 15 Figure 7: Cross-section showing ages of the geologic units in the region. The cross- section is an interpretation from seismic data conducted for oil exploration. The cross- section is a short segment from a longer seismic profile, marked in purple in the shaded relief, however, the exact location and boundaries of the cross section are not clear in the published literature. Something to note is the changing geology and faulting from north to south and the relative elevation of the basement. Adapted from Thomas et al. (1999). 16 The crust has an overall thickness of about 42 km in our area with a subcrustal Vp=8.20 km/s (Vs=4.7 km/s) under the Moho (Dimakov, 1968) (Figure 4). This profile interprets the Moho as a distinct boundary, but other studies (Abetov et al., 2015; Babadjanov et al., 1989, 1991; Sheikh-Zade, 1996) suggest a gradational Moho rather than a sharp one. The structure of this unit may be more complicated than what Dimakov (1968) shows in the profile. Mstislavskiy et al. (1980) argues that the Permian-Triassic Complex is highly deformed, folded, and uplifted in the northern part of Central Mangyshlak, but not elsewhere (Figure 8). This indicates that the seismic velocities of the Permian-Triassic Complex could change laterally due to the variable metamorphism in the uplift. The seismic velocities of the metamorphic parts of the Permian-Triassic Complex in this sub-region might be similar to those of the upper basement, thus might be misinterpreted. Further south, the structure changes again. Figure 4 show how the thickness of each unit changes as well the velocity (Dimakov, 1968). At the surface, South Mangyshlak begins with a 4 km thick Jurassic to Quaternary sedimentary basin with layer Vp=3.05 km/s (Vs=1.79 km/s). The underlying Permian-Triassic Complex also has a layer Vp=5.10 km/s (Vs=3.00 km/s) and a thickness of 3 km, thinning southward. The upper basement begins at the depth of 7 km and extends to 22 km with a total thickness of 15 km and a layer Vp=6.00 km/s (Vs=3.53 km/s). The lower basement has a total thickness of about 17 km and a layer Vp=6.80 km/s (Vs=4.00 km/s) ending at a depth of 40 km. The crust has an overall thickness of about 40 km in our area with a subcrustal layer Vp=8.00 km/s (Vs=4.7 km/s) under the Moho and the crust thins southward. Crustal structure of South Mangyshlak 17 Figure 8: Cross section showing folding in the western section of Central Mangyshlak from seismic imaging. Cross section from Mstislavskiy et al. (1980) where there is folding present in the Permian-Triassic Complex in Central Mangyshlak. The location cross section is shown in purple in the map. Cross section and legend are modified from Mstislavskiy et al. (1980). 18 DATA COLLECTION AND METHODOLOGY Seismic station deployment and data collection Ten broadband seismic stations were deployed throughout the Mangystau region in August – September 2016 (Mackey et al., 2017). All of the sensors were broadband Guralp CMG-3Ts paired with Reftek RT-130 dataloggers, borrowed from the IRIS/PASSCAL instrument center. All stations acquired data at 200 samples per second. Table 1 shows the coordinates of all the stations. Stations WKZ01, WKZ02, WKZ03, and WKZ04 were installed on metamorphosed sandstone/siltstone in Central Mangyshlak (Figure 1). Stations WKZ01 and WKZ02 were in the southeast while WKZ03 and WKZ04 were in the northwest. Stations WKZ05, 06, and 07 were deployed in South Mangyshlak, which is relatively flat and undisturbed. WKZ05 was installed on an unconsolidated layer of soil, while WKZ06 and WKZ07 were installed on limestone bedrock. WKZ05 is the westernmost station, while WKZ07 is the easternmost, near the border with Uzbekistan, with WKZ06 sited roughly midway between the two. Unfortunately, WKZ08 was erroneously removed by the Kazakhstan Border Patrol within a couple of days of installation. For this reason, little data are available so WKZ08 is not used in this study. WKZ09 was installed in a layer of unconsolidated sediment while WKZ10 was located on a limestone quarry outcrop. The deployment of the seismic stations used for this study was initially intended to determine the background levels of natural and anthropogenic seismic noise and determine if the region is suitable for the installation of permanent stations. These seismic Stations WKZ08, WKZ09, and WKZ10 were installed in North Ustyurt. 19 Table 1: Seismic stations used in this study. A list showing the location of each seismic station and the region where they were installed. 20 stations on the network were installed for about 20 days each, but not all at the same time (they all operated simultaneously for about 7 days). The first seismic station was installed on August 19, 2016 and the last one was demobilized on September 12, 2016. Once the data were collected, it was sent to the IRIS Data Management Center (DMC) for storage and cataloging, archived under network code “ZV”. The data used for this study was retrieved from the IRIS DMC and matched with 155 teleseismic events. Nine events were selected for use in this receiver function study (Table 2). Once the receiver functions were computed, they were inverted jointly with surface wave dispersion curves. The surface wave dispersion curves for each station location were derived from Pasyanos (2005) and provided by him. The method used to resolve for crustal structure was a joint inversion of receiver functions and surface wave dispersion curves from Julià et al. (2000). For this method, three inputs are needed for the inversion: receiver functions, surface wave dispersion curves for each station (or group), and an initial “a priori” velocity model. Receiver functions are 3-component, time series seismogram-like waveforms that show P to S converted phases caused by Earth structure (Langston, 1979). As the waves of an earthquake travel though the Earth, they will encounter velocity contrasts and create multiple conversions (P to S, and S to P waves), reflections, and refractions. From all the conversions produced, the P to S phases are useful since they are pulse-like arrivals and can be directly investigated on the radial component of a seismogram (Langston, 1979). We Receiver functions Methodologies 21 Table 2: Earthquakes used for receiver function estimations. The table lists earthquake parameters and the stations that produced an acceptable receiver function used for this study. The largest earthquake listed here produced a usable receiver function from every station. Refer to Figure 1 for a map of seismic stations and Table 1 for coordinates. 22 can calculate the thickness of the crust by measuring the time difference of each of these arrivals and using a previously known crustal velocity. Furthermore, the Earth structure under the seismic stations might be composed of different lithologies and have velocity contrasts that convert, reflect, or refract the energy into either P or S waves (Pp, Ps, PpPs, and PpSs + PsPs) that will create distinct peaks in the receiver functions. Figure 9 shows the different arrivals and their relative times with different color peaks for an idealized receiver function, and their respective travel paths inside the crust. For this method, we use teleseismic earthquakes that have occurred from 30 to 90° from the seismic station, since we want the angle of incidence to be as close as possible to vertical. Seismic events that occurred <30° from the station come in at a shallower angle and will travel longer paths within the crust to reach the seismic station, sampling different parts of the crust. As most seismic stations are deployed for long periods of time (many months or years), a great amount of data is usually available. In the case of this study, however, the seismic network was deployed for a maximum of ~20 days. To increase the available data, events with magnitudes as low as Mw 5.3 were used. Events in the data were correlated with National Earthquake Information Center catalogs (USGS NEIC) and were processed for the receiver functions. The first P arrivals were picked for all the earthquakes from the raw data, on all 3 components (vertical, north- south, and east-west). Once the arrivals were picked, the seismograms were cut 10 seconds before and 110 after the P arrival to eliminate unused data that might induce noise. The preparation of the data also involved 23 Processing Figure 9: Idealized ray diagram and receiver function of a teleseismic earthquake. a) A simplified ray diagram showing the converted phases used for the receiver function calculation. The solid lines are the converted P-waves while the dashed lines are the converted S-waves. b) The corresponding ideal radial receiver function of the model in (a). The color coding shows the relative arrival times of each converted phase. 24 Time Domain Iterative Deconvolution removing trends and tapering the signal as well as applying a bandpass filter of 0.05 to 8 Hz to eliminate high and low frequency data that represent noise. After the seismograms were prepared, they were rotated from the vertical, north-south, and east-west components to the radial and transverse components of motion. The radial component is the horizontal motion along the great-circle path from the earthquake to the receiver, and the transverse component of motion is orthogonal to the radial component in the horizontal plane. After the rotation, a time-domain iterative deconvolution approach (Ligorría and Ammon, 1999) was applied to eliminate the unwanted signal from the instrument response. Essentially, the results of the receiver functions are obtained by removing the instrument response from the radial component seismogram. The time-domain iterative deconvolution approach essentially creates a function that best models the Earth response using the data. First, the deconvolution strips out the largest peak from the observed data, followed by less important peaks. (Ligorría and Ammon, 1999). Ligorría and Ammon (1999) explain the process of deconvolving the waveforms iteratively by reducing the difference of the observed radial seismogram from a predicted signal generated by convolving an iteratively updated spike train (“earth response”) with the observed vertical component seismogram, using a least squares method. When the earth response and the vertical component are convolved, the radial component seismogram is produced and compared to the observed data. The first step for this deconvolution process is the cross-correlation of the vertical and the radial components to estimate the lag time of the first (and largest) peak in the spike-train. Once this spike train has the first peak, it is convolved with the vertical 25 component, and then subtracted from the radial component seismogram to observe how similar the signals are. The first iteration obviously won’t look like the radial component, so the process is repeated by cross-correlating the new estimate with the vertical component to find the next peak. The process is repeated until the difference between the processed signal and the observed radial component is minimized. For normal passive studies, the fit allowed for these seismograms is 85%, but due to our limited data, we allowed a fit of 75%. The iterative deconvolution program designed by Ligorría and Ammon (1999) allows as many iterations as desired, though were limited to 200 in this study After 200 iterations, the receiver functions that produced a misfit of 25% or less were further analyzed. The receiver functions were processed with 3 Gaussian filters (1.0, 2.5, and 5.0, which are similar to low pass filters of 0.56, 1.41, and 5.66 Hz respectively) to analyze different bandwidths. Most receiver function analyses use a Gaussian factor of 2.5 (Ligorría and Ammon, 1999), which provides more details in the signal (i.e. shows more conversions from structures within the crust), but in this study only a few receiver functions had a fit better than 75% with the higher Gaussian filters. The higher Gaussian filters (2.5 and 5.0) show more details about the Earth structure but can also introduce more noise due to the instability of the method (Figure 10). Therefore, all of the receiver functions shown in this study were processed using a Gaussian filter of 1.0, where only low frequency data are used. Stacking Ideally, one wants the earthquakes to occur from all different azimuths in order to reduce any bias in signals coming only from limited direction. This problem persists in earthquake seismology since earthquakes are not evenly distributed. Some directions 26 09/09-11:53:15 Gaussian = 5.0 Gaussian = 2.5 Gaussian = 1.0 Pp Ps PpPs 5 5 PsSs + PsPs 15 20 25 0 10 Figure 10: Example of receiver functions with different Gaussian filters from earthquake detected in station WKZ10. The 3 receiver functions are from the same earthquake (earthquake information on box on top-right corner) and showing the ideal phases for a receiver function in the waveform processed with a Gaussian filter of 1.0. As the value of the Gaussian filter increases, more information appears in the waveform or more noise appears. The Pp peak splits up into multiple peaks in the higher order filters. Time (s) 30 27 might show evidence of structure that others won’t detect, so the signals would show slightly different information. Due to this problem, and to increase the signal to noise ratio, we stacked the waveforms from each station to obtain an average of crustal thickness. For each stack, we could estimate crustal thickness by multiplying an assumed crustal velocity times the time difference between the Pp and Ps. The two main problems here are: 1) The crust is not a single layer with a single velocity, though this is the assumption that is made, and 2) the solution is non-unique (Ammon et al., 1990), where we have a trade-off between velocity and thickness that can give us the same time calculated by the receiver function (Ps - Pp arrival times). For this reason, we decided to approach the problem of resolving the crustal structure within the crust, as well as the overall thickness, using the method of Julià et al. (2000), which is a joint inversion of receiver functions with surface wave dispersion observations. Similar to receiver functions, the dispersion of surface waves is also sensitive to shear velocity. The longer period surface waves travel faster and penetrate deeper than the shorter period waves, and both are affected by the shear velocity in the crust. Unlike receiver functions, surface wave dispersion constrains the averages of shear wave velocity at which these periods travel. The surface wave dispersion curves for each individual station used in this study were provided to us from Pasyanos (2005) from a variable resolution model, where the study can resolve details up to 1°. The surface wave dispersion measurements were group velocities for different periods (from 10 s to 100 s in intervals of 5 s) of both Love and Rayleigh waves. Surface Wave dispersion 28 The dispersion of Rayleigh waves is calculated differently than Love waves due to their distinct particle motions. A more thorough explanation for obtaining the dispersion curves and surface wave tomography process is given by Pasyanos et al. (2001). Normally, a narrow-band Gaussian filter is applied to the broadband vertical component over many periods to obtain the Rayleigh wave dispersion, which is subsequently applied to the rotated transverse component to find the Love wave dispersion. For each one of these observable periods, the maximum amplitude in the filtered waveform is picked, and the arrival time is used to compute the group velocity since it is the most representative peak of that group. The group velocity is then matched to a period using the PREM model from Dziewonski and Anderson (1981), which shows the relationship of period vs group velocity. After finding the periods of all the group velocities calculated, the surface wave dispersion is used to estimate lateral group velocity variations by doing a travel-time inversion with a conjugate gradient method. The relationship between Rayleigh/Love dispersion and Earth structure can be obtained by calculating sensitivity kernels (Figure 11), which show the sensitivity of different period waves to different depths within the crust. These are obtained by calculating the sensitivity of the seismograms to any perturbations in the parameters in the earth models, using a crustal velocity model (Randall, 1994). To obtain a reliable surface wave dispersion model, there is a minimum threshold of observations that should be included in the period calculations. Pasyanos (2005) points out that periods having less than 5,000 picked arrivals are not reliable, which is the case for period of 5 s in his model. The 10 s period had less than 10,000 picked arrivals, compared 29 Figure 11: Shear wave sensitivity kernels for Love and Rayleigh waves. The sensitivity kernels from the Pasyanos (2005) model are calculated with the velocity model on the left of the figure. The velocity model shows typical shear (dashed line) and compressional (solid line) velocities of continental crust. Love wave sensitivity kernels are shown center and Rayleigh wave sensitivity kernels at right. Figure from Pasyanos (2005). 30 to about 23,000 for the periods with most picks, and was also eliminated from our inversion process. The Love waves have much lower sensitivity than Rayleigh waves and they don’t sample as deep. The longer period Rayleigh waves are more sensitive to deeper structure, as opposed to Love waves, where all of the periods up to 60 s are sensitive to shallow structure. The depth resolution for each type of wave is different as a result of the physics of particle motions. The periods we focus on are the 15 – 60 s Love waves, which are much more sensitive to upper structures and basins. Even though the Love waves have a lower signal to noise ratio, we have many periods that sample very similar depths and we can thus obtain a reliable image of the upper to middle crust. To look at deeper crust, we focus on the 15 – 20 s Rayleigh waves to observe deep basins, and 30 – 50 s for crustal thickness (Moho depth). Pasyanos (2005) shows that 60 s Rayleigh wave group velocities can also see crustal velocities present in the Mangystau Region. The surface wave dispersion does not change significantly from station to station in the Mangystau Region (Figure 12). The Rayleigh wave dispersion curves for all stations seem to be constant throughout, but the one for WKZ07 is faster than the rest. There seems to be no significant difference between the stations that were installed in different lithologies. On the other hand, the group velocities of the Love wave dispersion curves differ more for lower periods. All the curves follow the same pattern, but the 20 s period waves from WKZ07 are about 0.29 km/s faster than the next faster station, WKZ01, and the difference becomes smaller with increasing periods. The fastest station is WKZ07 for both Rayleigh and Love waves in the periods of interest. There is no significant difference in patterns of curves or velocities that gives us reason to prefer using the curve of one station 31 Figure 12: Surface wave dispersion curves for all the seismic stations. Graph showing the Rayleigh (top) and Love wave (bottom) dispersion curves with a logarithmic scale for the Period axis. The two types of waves follow different patterns since they are affected differently at different depths, but they are mostly consistent throughout the stations. 32 vs the other. Once we have receiver functions and surface wave dispersion datasets, we can apply the method which combines them and inverts for a surface wave velocity model. This can be done since both datasets contain similar information calculated with different parameters. Initial shear wave velocity structure models The third input needed for the joint inversion is the initial shear wave, or “a priori”, velocity model. We used velocities from Dimakov (1968), Mstislavskiy et al. (1980), and Sheikh-Zade (1996) to construct these initial models, because they represent measurements obtained through other seismic methods. The most relevant study was Dimakov (1968) because he obtained measurements of the same area we are studying. We created 3-layer crustal velocity models with the following parameters: upper layer with Vs ranging 1.7 – 3.30 km/s, middle layer with Vs ranging 2.8 – 3.6 km/s, and lower layer with Vs ranging 3.3 – 4.00 km/s. The upper layer represents the Jurassic to Quaternary and Permian-Triassic sedimentary units in the shallow crust, the middle layer represents the granitic-metamorphic unit of the basement, and the lower layer represents the basaltic unit of the basement from Dimakov (1968). We then shifted the boundaries between each layer to make them thinner or thicker relative to each other. Finally, we shifted the Moho from 30 to 42.5 km to see how the final models would change with different crustal thicknesses. Figure 13 shows a representation of how the layers were shifted and what ranges of velocities we used. Once we have our 3 main inputs, we used a linearized joint inversion of receiver functions and surface wave dispersion to obtain information on geology of the upper Joint inversions of receiver functions and surface wave dispersion 33 Figure 13: Diagram showing the range of “a priori” shear wave velocity structure models used. The velocity model is used as one of the main inputs for the joint inversion. Various velocity models were generated by changing the layer thicknesses, velocities, and Moho depth. The values are next to the name of each section. The red values to the left side of the velocity line are Vp values (in km/s) corresponding to those Vs values by using a Vp/Vs ratio of 1.7. 34 mantle and crust (Randall, 1994; Julià et al., 2000; Maceira et al., 2009; Chai et al., 2015). Receiver functions and surface wave dispersion curves are inverted together because receiver functions constrain information on shear velocity contrasts in the crust and their depths (Ammon et al., 1990), while surface wave dispersion constrains the averages of the absolute shear velocity at different depths (Julià et al., 2000). We can benefit from this method since the calculations of both datasets are independent from each other and the surface wave dispersion model adds a second constraint on the inversion. Inverting both methods together will constrain the shear velocity structure better than either method individually. However, these methods represent different physical parameters and have a different sets of data points (receiver functions are time vs amplitude and surface wave dispersion are period vs group velocity), but the inversion algorithm compensates for that by equalizing both datasets and giving different weights to each as desired. The equalization is achieved by dividing the desired weight of each dataset by the number of data points and the variance through a differential damped least squares calculation. Typical values for variances of receiver functions are 0.01 s-1 and for dispersion measurements these values are 0.05 km/s according to Julià et al. (2000). The influence of the initial “a priori” velocity model in the inversion process can be manipulated. We can decide how much weight the initial velocity model has on the inversion process, where a higher a priori weight will have more influence in the inversion. The values given to this influence are arbitrary, and the best way to manage it is by giving many values to it and observe the outcome. For this study, we tested four different values [0.3, 0.5, 0.7, and 1.0] for the influence of the a priori velocity model in the inversion process. The lowest value represents the least influence of the “a priori” model in the 35 inversion, resulting in a final model mostly computed from only the dispersion and receiver function inversion. The highest value though, produces a final model that looks almost exactly like the a priori velocity model. Second, we can control the relative weight of each dataset used in the final result. We can assign the weight of each dataset by specifying a value ranging from 0 to 1, where 0 is 100% weight given to receiver functions and 1 is 100% weight given to surface wave dispersion. We use values of 0.25 (75 % receiver functions, 25% dispersion) and 0.50 (50% each dataset) for specific cases, which are discussed in the results. The stations deployed in Central Mangyshlak sampled a different structure (uplifted basement) than the rest of the stations; however, the dispersion model of Pasyanos (2005) does not show any evidence of such an uplift, only basin velocities. The intermediate and long period waves can’t resolve structures at the scale of Central Mangyshlak, and the short period waves that would, are not reliable due to lack of data from the dispersion model (this is also the case for all the station locations). Therefore, the datasets of stations WKZ01 – WKZ04 were weighted with both values (0.25 and 0.50). This is important because we can give preference to one dataset over the other if one of the datasets lacks, or has erroneous, data points and can be compensated by the other. Finally, we smooth the final velocity model to eliminate any sharp boundaries that can be artifacts of noise in the datasets. This value is also arbitrary, but its effect can be determined by trying different numbers in the inversion process. The five values tested were 0.15, 0.7, 1.50, 3.00, and 5.00. The first value did not give any significant smoothness, and the highest value didn’t converge to a final velocity model in most cases. The most reliable values, and the main ones used in this study, are 0.70 and 1.50, which smooth the 36 final result, but still show some smaller structure within the crust (Figure 14). Once has to be careful increasing the smoothing, because it makes sharp boundaries look smooth, like the Moho. Once the inversion process is complete, a new receiver function and new surface wave dispersion curve are computed from the final velocity model, so we can observe the misfit between the observed and final results. We adopt the values of variance for a regular receiver function and surface wave dispersion to be the allowable threshold for misfit between the observed and final result of the joint inversion. Typically, other studies allow a misfit of 0.05 between observed and final dispersion curves, and 0.01 between observed and final receiver functions (Juliá et al., 2000). This misfit is found by calculating the sum of the difference between the samples of the observed data and the final (processed) data. However, we allowed the misfit to increase as long as the misfit from the other dataset is closer to the allowed misfit value, and the final velocity model makes sense for the following reasons: the periods that we are using don’t sample the scale in the upper crust that we want, and second, the data are limited meaning that we don’t have a very good data coverage of the region. Finally, the final velocity model is used for interpretation if the structure and velocities are somewhat similar to those from previous studies, and the misfit of the receiver functions and surface wave dispersion curves are near the threshold. To summarize, the joint inversion is computed by equalization of both datasets. Since each dataset has its errors (or instability in the models) we smooth the final velocity model resulting in an overall velocity structure with gradients instead of discrete layers. The inversion depends on an initial velocity model, and each iteration of the inversion updates the final velocity model and creates a receiver function and dispersion curve based 37 Smoothness = 0.15 Smoothness = 0.70 Smoothness = 1.50 Smoothness = 5.00 Smoothness = 3.00 Figure 14: Comparison of different smoothness factors from 0.15 to 5.00. The 5 models show an initial and final velocity models for WKZ01. Each final velocity model is different based on the smoothness constraint applied to that iteration (smoothness value below each corresponding graph). The smoothness of 0.15 shows many layers that are artifacts and this jagged pattern decreases with increasing smoothness. The smoothness of 5.00 deviates the final velocity model completely, hence it was not used in the processing. The best results come from a smoothness of 0.70 and 1.50. This is not a definitive final result, but rather to show a comparison between the impact of smoothness constraints. 38 on current velocity structure, until the 8th iteration. Julià et al. (2000) point out that this number of iterations is enough for the program to converge on a result. If we run more than 8 iterations, the results will become better, but the improvement will be insignificant. More advanced methods have been used in the last decade but need much more data and need additional types of data. A Monte Carlo approach can be modeled to infer the crustal structure (Shen et al., 2013), but we need more data in order to achieve this. Another method would be jointly inverting receiver functions and surface wave dispersion with gravity data (Maceira et al., 2009; Chai et al., 2015), but we don’t have gravity data. The code we used was the same used to invert the seismic data with gravity from Chai et al. (2015), but we omitted the option to invert with gravity and only used the joint inversion of receiver functions with surface wave dispersion. 39 RESULTS FOR RECEIVER FUNCTIONS AND LINEARIZED JOINT INVERSION OF Receiver functions results RECEIVER FUNCTIONS AND SURFACE WAVE DISPERSION In this section, we present the receiver functions computed for all of the stations in this seismic experiment, except for WKZ08. Receiver functions were computed for each individual station though the number of receiver functions with a fit higher than 75% varied. It is important to mention that lateral homogeneity is assumed for every station. We only used receiver functions calculated having a Gaussian filter factor of 1.0 because higher filters introduced noise to the waveform. From 151 teleseismic events, 9 produced receiver functions having a misfit of less than 25% (fit of 75%). Figure 15 compares all of the receiver functions produced in the process. We can see the Pp and Ps arrivals clearly, but there is no obvious trend for the PpPs or PsSs+PsPs phases. One significant observation is the width of the Pp arrivals; the arrivals for stations WKZ01 through WKZ04 are narrower (~ 3 seconds) than the ones for the rest of the stations (~ 5 seconds). This is likely due to the fact that stations WKZ01 – 04 were placed in the Central Mangyshlak uplift, where there is an exposed basement layer of metamorphosed sediments. The rest of the stations were installed in North Ustyurt and Southern Mangyshlak. The rest of the stations have a wider Pp arrival since they include a conversion caused by the base of the sedimentary basins that can’t be resolved with the Gaussian filter of 1.0. The Ps arrivals don’t line up perfectly but can be followed throughout most of the waveforms. One can make interpretations about the crustal thickness using only receiver functions, but as Ammon et al. (1990) explain, the solution is non-unique due to the trade- 40 Figure 15: Individual receiver functions from all the stations. The horizontal axis shows the time of the function and the vertical shows the relative amplitudes to compare the receiver functions. The station name for each receiver function is shown on the left side of the figure. Peaks and troughs are filled to observe the phases easier. The Pp arrival is narrower in stations WKZ01 – 04 and is wider in the rest of the stations. The Ps phase can be followed through very well, but later reverberations can’t be correlated without stacking. 41 off between assumed velocity and depth. Since each station has a different number of results (WKZ04 has one versus WKZ05 with five receiver functions), and the shear wave velocity (Vs) of the region is not well constrained, we do not attempt to conclude a Moho depth using just receiver functions. Considering the limited number of receiver functions available and the difficulty in identifying the PpPs and other phases on individual traces, we stacked the receiver functions by station and/or groups of stations. With stacked receiver functions, the signal to noise ratio increases and low-amplitude phases become clearer. Beginning with North Ustyurt, four receiver functions were produced for station WKZ09 where the Pp, Ps, PpPs, and PsSs+PsPs phases are very clear (Figure 16). The PpPs phase in the “09/09 – 11:53:15” receiver function arrives earlier than those of the other receiver functions, which causes that phase in the stack to be broader and have a double peak. Four receiver functions were computed from station WKZ10 (Figure 17). The second and third receiver functions (top to bottom) in Figure 17, show the Pp phase skewing, as if it was being smoothed with another phase. The first receiver function shows PpPs and PsSs+PsPs phases that arrive later than in other receiver functions from the same station. The stack for WKZ10 has very clear Pp and Ps, but unclear PpPs and PsSs+PsPs phases. The PpPs phase has a double peak and covers the PsSs+PsPs phase. Since both stacks contain some noise that comes from the shift in PpPs and PsSs+PsPS phases, we create a group stack of stations WKZ09 and WKZ10 (Figure 18); they are in the same tectonic unit, and they have very similar receiver functions. If we compare 42 Pp Ps PpPs 09/09-11:53:15 09/05-22:54:03 08/29-04:29:57 08/24-10:34:54 PsSs + PsPs 5 0 5 10 15 20 25 30 Time (s) Pp Ps PpPs PsSs + PsPs 5 10 5 0 15 Figure 16: Individual and stacked receiver functions of station WKZ09. (top) Four receiver functions computed from station WKZ09 showing the clear Pp and Ps arrivals and reverberations. (bottom) The stack created with the 4 receiver functions from the figure on top. The shift in PpPs causes a double peak in the stack. The events used to compute receiver functions are shown in order in the box in the top-right corner of the upper graph. Time (s) 30 25 20 43 Pp Ps PpPs 09/09-17:17:54 09/09-11:53:15 09/05-22:54:03 08/29-04:29:57 PsSs + PsPs 5 0 5 10 15 20 25 30 Time (s) Pp Ps PpPs PsSs + PsPs 15 10 Figure 17: Individual and stacked receiver functions of station WKZ10. (top) Four receiver functions computed from station WKZ10 showing all the phases of an ideal receiver function with some noise. (bottom) The stack created with the 4 receiver functions from the top image showing a higher signal to noise ratio. The events used to compute receiver functions are shown in order in the box in the top-right corner of the upper graph. Time (s) 20 25 30 5 0 5 44 Figure 18: Stack from group of stations WKZ09 and WKZ10. Graph showing the stack created by the receiver functions of stations WKZ09 and WKZ10, from figures 16 and 17, respectively. Note the cancelling of noise after the Ps phase, where the individual stacks had either a broad peak (WKZ09) or a discrete double peak in the PpPs phase (WKZ10) (Figure 17). 45 their stacks, there is much overlap and the minor peaks from each stack cancel each other. This leaves us with a clean receiver function from the grouping. In Central Mangyshlak, three receiver functions were produced for WKZ01 (Figure 19) with a consistent Pp width and a Ps phase at 5 seconds. The PpPs and PsSs+PsPs phases in these waveforms arrive between16 and 18 seconds and are very clear. The stack of the three receiver functions produces a signal very similar to the individual original receiver functions, suggesting consistency and repeatability of the receiver functions. Three receiver functions were produced for WKZ02 (Figure 20) with a less consistent Pp width. The Ps arrivals are consistent with those of WKZ01 and with each other, except for the “09/05– 22:54:03” receiver function, which reduces the amplitude of the Ps arrival in the stack. The stack from this station shows less noise than the individual receiver functions. The arrival times of all phases are similar, but the amplitudes changed because they do not overlap throughout the waveforms. Since WKZ01 and WKZ02 have only three receiver functions each, are in the same tectonic unit, are only 50 km apart, and the Pp and Ps arrival times are less than 0.5 seconds apart, we grouped these two stations together and made a single stack of their receiver functions (Figure 21). Both stations show very similar Pp, Ps, PpPs, and PsSs+PsPs arrival times. The Pp phase in the receiver functions of this group have roughly the same width, but the amplitudes change between both stacks. The stack of WKZ02 might have decreased in amplitude because of the waveform that had the missing Ps phase. Four receiver functions were computed for station WKZ03 (Figure 22). The Ps phase is very clear in the “09/10– 12:27:33” and “08/29– 04:29:57” receiver functions, but its amplitude is significantly smaller in the other two receiver functions. However, the 46 Pp Ps 08/29-04:29:57 08/24-10:34:54 08/24-01:36:32 PpPs PsSs + PsPs 5 0 5 10 15 20 25 30 Time (s) Pp Ps PpPs PsSs + PsPs 5 0 15 Time (s) Figure 19: Individual and stacked receiver functions of station WKZ01. (top) Three receiver functions computed from station WKZ01 showing all the phases of an ideal receiver function with some noise. (bottom) The stack created with the 3 receiver functions from the top image showing a higher signal to noise ratio. The events used to compute receiver functions are shown in order in the box in the top-right corner of the upper graph. 30 20 25 5 10 47 Pp Ps PpPs 09/05-22:54:03 08/29-04:29:57 08/24-10:34:54 PsSs + PsPs 5 0 5 10 15 20 25 30 Time (s) Pp Ps PpPs PsSs + PsPs 5 0 5 10 15 Time (s) Figure 20: Individual and stacked receiver functions of station WKZ02. (top) Three receiver functions computed from station WKZ02 showing all the phases of an ideal receiver function with some noise. (bottom) The stack created with the 3 receiver functions from the top image showing a higher signal to noise ratio. The events used to compute receiver functions are shown in order in the box in the top-right corner of the upper graph in order. 25 20 30 48 Figure 21: Stack from group of stations WKZ01 and WKZ02. Graph showing the stack created by the receiver functions of stations WKZ01 and WKZ02, from figures 19 and 20, respectively. All the phases are clear in the stack since the noise was reduced. 49 Pp Ps 09/10-12:27:33 09/09-11:53:15 09/05-22:54:03 08/29-04:29:57 PpPs PsSs + PsPs 5 0 5 10 15 20 25 30 Time (s) Pp Ps PpPs PsSs + PsPs 5 0 20 25 30 5 15 10 Time (s) Figure 22: Individual and stacked receiver functions of station WKZ03. (top) Four receiver functions computed from station WKZ03 showing all the phases of an ideal receiver function with some noise. These waveforms don’t show the PpPs phase very well, but it is labeled where it should show up based on other stations. (bottom) The stack created with the 4 receiver functions from the top image showing a higher signal to noise ratio. The events used to compute receiver functions are shown in order in the box in the top-right corner of the upper graph. 50 arrival time and width of the Pp phase is consistent in all the receiver functions. The PpPs phase is not very clear in any of these receiver functions but the PsSs+PsPs phase arriving at about 20 seconds is constant in all the waveforms. The stack shows the PpPs and PsSs+PsPs phases within the time window of 16 – 20 seconds. Although WKZ04 was an exceptionally quiet station, it only produced one receiver function (Figure 23). This receiver function clearly shows all the phases that are in an ideal receiver function in the same time window as WKZ03. However, there is lower amplitude noise in the time between the Ps and PpPs phase arrivals. For the same reason we grouped stations WKZ01 and WKZ02, we grouped stations WKZ03 and WKZ04. The distance between these stations was only 60 km. The stack created by averaging the stacks of WKZ03 and WKZ04 looks much more like the one of WKZ03. (Figure 24). However, there is still a lot of noise in between the Ps phase and PpPs phases. In South Mangyshlak, five receiver functions were produced for station WKZ05 (Figure 25). The Pp and Ps arrivals are constant throughout all of the receiver functions as well as the PsSs+PsPs phase with its arrival time between 15 to 20 seconds. The stack of WKZ05 looks almost the same as its individual receiver functions (Figure 25, bottom graph). Four receiver functions were produced for station WKZ06 (Figure 26). An extra peak between the Ps and PpPs phases can be traced throughout all the receiver functions of WKZ06 and is higher than the Ps arrival in some cases. This extra peak might be indicative of another major velocity contrast within the crust. After stacking the receiver functions, 51 Pp Ps PpPs PsSs + PsPs 08/29-04:29:57 0 5 5 10 Figure 23: Receiver function of station WKZ04. This plot shows the only receiver functions computed from station WKZ04. This station can’t be used individually for further processing since there are not enough waveforms to find a more unique solution for crustal thickness. The event used to compute the receiver function is shown in the box in the top-right corner of the graph. Time (s) 30 15 20 25 52 Figure 24: Stack from group of stations WKZ03 and WKZ04. Graph showing the stack created by the receiver functions of stations WKZ03 and WKZ04, from figures 22 and 23, respectively. The Pp and Ps arrivals are clear, but he signal has a negative amplitude between phases Ps and PpPs where it is difficult to map were the PpPs and PsSs+PsPs phases are. 53 Pp Ps PpPs PsSs + PsPs 09/05-22:54:03 08/29-04:29:57 08/24-13:48:45 08/24-10:34:54 08/20-09:01:26 5 0 5 10 15 20 25 30 Time (s) Pp Ps PpPs PsSs + PsPs 5 0 5 10 15 Time (s) Figure 25: Individual and stacked receiver functions of station WKZ05. (top) Five receiver functions computed from station WKZ05 showing all the phases of an ideal receiver function with some noise. (bottom) The stack created with the 5 receiver functions from the top image showing a higher signal to noise ratio. The events used to compute receiver functions are shown in order in the box in the top-right corner of the upper graph. 30 20 25 54 Pp Ps PpPs PsSs + PsPs 09/05-22:54:03 08/29-04:29:57 08/24-13:48:45 08/24-10:34:54 5 0 5 10 15 20 25 30 Time (s) Pp Ps PpPs PsSs + PsPs 5 0 5 10 15 Figure 26: Individual and stacked receiver functions of station WKZ06. (top) Four receiver functions computed from station WKZ06 showing all the phases of an ideal receiver function with some noise. (bottom) The stack created with the 4 receiver functions from the top image showing a higher signal to noise ratio. The events used to compute receiver functions are shown in order in the box in the top-right corner of the upper graph. Time (s) 25 20 30 55 the Ps has a higher amplitude than the extra peak, which should be the case due to the presumed velocity contrast from the Moho being higher than anything within the crust. Since the extra peak between the Ps and PpPs phases appears in three of the receiver functions we speculate the extra peak might be caused by structure. This arrival does not appear on WKZ04, which can imply that it was caused by local structure and is not directional dependent. More data and further studies are needed in order to conclude on the source of this peak. WKZ05 and WKZ06 are further away from each other than for the other groups produced (120 km), but they are also in the same tectonic unit and the phase arrival times are similar to each other. The stack created by the receiver functions of WKZ06 and WKZ05 increases the signal to noise ratio significantly and the phases are more evident than either individual stack (Figure 27). This group has the most receiver functions (total of 9) compared with the other groups. Finally, station WKZ07 is isolated in South Mangyshlak. The nearest station WKZ01 is 140 km away but doesn’t share the same tectonic setting. WKZ06 shares the same tectonic setting but is 180 km away. Therefore, WKZ07 was not grouped with other stations. The receiver functions from this station are quiet, with minimal noise, the Pp and Ps phases are clear and match each other, but the PpPs and PsSs+PsPs phases can be traced in only one receiver function (Figure 28). The “08/29 - 04:29:57” receiver function almost looks like an ideal receiver function, but the Ps phase arrives earlier in comparison to the other receiver functions. 56 Figure 27: Stack from group of stations WKZ05 and WKZ06. Graph showing the stack created by the receiver functions of stations WKZ05 and WKZ06, from figures 25 and 26, respectively. Pp and Ps phases are clear but PpPs and PsSs+PsPs are more difficult to interpret as they have similar amplitude than the noise. 57 Pp Ps 08/29-04:29:57 08/24-10:34:54 08/24-01:36:32 PpPs PsSs + PsPs 5 0 5 10 15 20 25 30 Time (s) Pp Ps PpPs PsSs + PsPs 5 5 15 20 25 0 10 Time (s) Figure 28: Individual and stacked receiver functions of station WKZ07. (top) Three receiver functions computed from station WKZ07 showing all the phases of an ideal receiver function with some noise. Only the top receiver function shows the PpPs and PsSs+PsPs phases and the label of these phases shows where they should be on the other waveforms as well. (bottom) The stack created with the 3 receiver functions from the top image showing a not so obvious PpPs and other reverberated phases. The events used to compute receiver functions are shown in order in the box in the top-right corner of the upper graph in order. 30 58 Joint inversion of receiver functions with surface wave dispersion results Once the receiver functions are prepared, stacked, and grouped, they were used for the joint inversion process. Each station’s stack was inverted individually, as well as the groups of stations. The joint inversion of receiver functions (RF) and surface wave dispersion curves were computed by testing different parameters. Since most of the individual stations did not have enough data to determine a final velocity model, we used the groups of stations to find the shear wave velocity structure of each tectonic unit. For the inversion process, the following parameters were controlled and manipulated: the weight of the “a priori” velocity model, the initial velocity models, the influence of each dataset (RF and dispersion), and the smoothing factor. Some groups had special circumstances where the parameters had to be manipulated differently based on their noise and initial results and are explained below. We processed the joint inversion of RF and surface wave dispersion with 48 a priori models to simulate different structural scenarios (see page 33 for details). The Moho depths used varied between 35 to 42.5 km in increments of 2.5 km, where the velocities changed from lower crustal to upper mantle velocities. All of the receiver functions were inverted using all of the a priori models. In theory, all a priori models used should converge to the same model. However, not all of the a priori models converged and some that did, yielded inconsistent results. To determine the relative weight of the receiver functions and dispersion curves used in the inversion, we used two different scenarios. The first scenario was to invert the groups in Central Mangyshlak using weights of 75% RF and 25% dispersion. The Central Mangyshlak stations WKZ01 through WKZ04 were installed in a structure smaller than 59 what the surface wave dispersion studies (Pasyanos, 2005) can resolve, thus the weght assigned to the influence of the dispersion model is low. The receiver functions produced from these stations contain some information about the shallow structure, as shown in Figure 15, where the Pp phases are much narrower than those of the of the stations installed over deep basins. The second scenario was to invert the other groups, not in Central Mangyshlak, with a 50% weight on RF and 50% on surface wave dispersion. This is because both receiver functions and surface wave dispersion results show evidence of sedimentary basins; receiver functions show a wider Pp phase (see page 26), and longer period surface waves can image them (based on sensitivity kernels). RF and surface wave dispersion curves are weighted equally as they are used with equal confidence. Finally, as mentioned earlier, we focus more on the models inverted with smoothness values of 0.5 and 0.7. The higher values smooth the final velocity model more and decrease the number of discrete layers (Figure 14). Moreover, this produces final velocity models with a velocity gradient instead of sharp discontinuities. The smoothing used for the final velocity models depends on how noisy the model came out to be. The depth of the boundaries in the crust were determined by the inflection points, or slope breaks, in the final velocity models. This is because the sharp steps in the final velocity model are being smoothened. The Moho boundary is defined at Vs=4.4 km/s, since this value is representative of upper mantle velocities. From more than 3,500 models that were processed, we only show the results with the lowest misfit between observed data and final modeled data (RF and surface wave dispersion curve calculated from the final velocity model), and/or if the final velocity model makes sense with the geology of the area. 60 The first group contains stations WKZ09 and WKZ10 located on North Ustyurt. Using an 0.5 a priori weight, a 50 – 50% weight from each receiver function and dispersion, and a smoothness of 1.50, the joint inversion of this group produces a good, but smooth, velocity model (Figure 29). This final velocity model begins with two upper sections starting with a velocity of 1.33 km/s at the surface and steadily increasing to 2.13 km/s at 3 km depth, followed by a section increasing from 2.64 km/s to 3.34 at 10.5 km depth. The following section does not change much, with velocity increasing from 3.51 km/s to 3.61 km/s velocity at 20.5 km depth. The lower section increases from 3.76 km/s to 3.88 km/s at 35.5 km depth, where it is followed by a gradient to 41 km depth with velocities increasing from 3.98 km/s to 4.67 km/s within the gradient. If the inversion was not smoothed, the final model would probably have a Moho step somewhere in the gradient from 36 to 41 km depth. The misfits from the final and the observed data are 0.0942 for Rayleigh waves and 0.1421 for Love waves. The Rayleigh waves with periods of 20 – 45 s from the final model curve fit those from the observed data as well as the final model 25 – 65 s period Love waves fit the observed data. The receiver function misfit is 0.0269 where the final RF Pp and Ps amplitudes are lower and the troughs between them are higher. In Central Mangyshlak, the next group of stations WKZ01 and WKZ02 are used (Figure 30). The inversion with a priori weight of 0.5, weight of 75% RF and 25% dispersion, and a smoothness of 0.70 produced the best result for this grouping. The 0.70 smoothness was used because the lower smoothness resulted in unrealistic low velocity layers in the lower crust, and the higher smoothing completely smoothed the model. 61 Figure 29: Joint inversion result for group of stations WKZ09 and WKZ10. The box at left shows the initial a priori velocity model (gray) compared to the final velocity model (red). The green dashed lines represent the slope breaks that mark the boundary of the sections. From the final velocity model, a final receiver function (bottom right) and surface wave dispersion curves (top right) were calculated. The “initial” results are calculated from the a priori model to show what the results would be if the model was true. The final velocity model shows two upper sections of 3 km and 7.5 km, a 11.5 km middle section, a 16 km lower section and a 5 km thick gradient at the bottom of the crust indicating an overall thickness of 36 – 41 km. 62 Figure 30: Joint inversion result for group of stations WKZ01 and WKZ02. The box at left shows the initial a priori velocity model (gray) compared to the final velocity model (red). The green dashed lines represent the slope breaks that mark the boundary of the sections. From the final velocity model, a final receiver function (bottom right) and surface wave dispersion curves (top right) were calculated. The “initial” results are calculated from the a priori model to show how the results would be if the model was true. This final velocity model shows two upper sections of 3 km and 7.5 km, a 12.5 km thick middle section, a 12.5 km thick lower section, followed by a 7.5 km thick gradient at the base of the crust with an overall thickness of 36 – 43 km. 63 The first section has a velocity of 1.97 km/s at the surface and steadily increases to 2.46 km/s at 3 km depth. Another upper crustal section from 3 km to 10.5 km depth is present with velocities increasing from 2.85 to 3.25 km/s. The middle section of the crust follows with velocity of 3.49 km/s that very steadily increases to 3.57 km/s at 23 km depth. The lower section is found from 23 to 35.5 km, with velocity increasing from 3.68 to 4.03 km/s. The base of the crust has a velocity gradient of 4.21 km/s at 35.5 km depth to 4.64 km/s at 43 km depth, where the Moho is thought to be. It is unclear if the Moho is gradational or if it is smoothened by the 0.70 parameter. The misfit between observed and final surface wave dispersion curves is 0.0324 for Rayleigh waves and 0.2664 for Love waves. The Rayleigh wave observed and final dispersion curves match at 20 to about 70 seconds, but the final velocities slightly overestimate the observed velocities in the intermediate periods. The Love wave dispersion curves match from periods of 25 seconds to 70 seconds. The misfit between the observed and final receiver functions was 0.0268. The amplitude of the Pp phase in the final receiver function is about 0.125 s-1 lower than the observed but the rest of the waveform matches almost perfectly up to 19 seconds. The inversion of the grouping of stations WKZ03 and WKZ04, also in Central Mangyshlak, show inconclusive results. The best model produced was inverted with a priori weight of 0.3, weighting of 75% for receiver functions, and 0.7 smoothness (Figure 31). In the upper section, the velocity begins at 0.97 km/s, increases to 2.31 km/s until 5 km, then decreases to 2.01km/s for the next 2 km, and suddenly increases to 2.7 km/s. Based on the geology of the region and previous studies, no one has found a layer of lower velocity material under sediments, nor have low velocities such as 0.97 km/s been 64 Figure 31: Joint inversion result for group of stations WKZ03 and WKZ04. The box at left shows the initial a priori velocity model (gray) compared to the final velocity model (red). From the final velocity model, a final receiver function (bottom right) and surface wave dispersion curves (top right) were calculated. The “initial” results are calculated from the a priori model to show how the results would be if the model was true. This final velocity model shows no convergence on any section and is not reliable. 65 reported. The surface wave dispersion misfits are way beyond the threshold at 0.3471 for Rayleigh waves and 1.0191 for Love waves. The observed receiver function does not match the final receiver function, except for the Ps phase, with a total misfit of 0.0619 between both waveforms. Therefore, the final velocity model produced from the group of stations WKZ03 and WKZ04 does not appear to be reliable. There were two models created for South Mangyshlak. We present the first model, from grouping of stations WKZ05 and WKZ06, which it was inverted by a priori weight of 0.5, a weight of 50 – 50% to receiver functions and dispersion, and a smoothness of 0.70 (Figure 32). Like other groups, there are two upper sections where the top starts at 1.35 km/s and increases to 2.63 km/s to 5.5 km depth, where the other section increases from 3.01 to 3.25 at 10.5 km depth. The middle section has velocity slowly increasing from 3.49 km/s to 3.58 km/s at 23 km depth, followed by the lower crustal section with velocity increasing from 3.67 km/s to 3.96 at 33 km depth. There is a velocity gradient at the bottom of the crust from 33 km to 40.5 km depth, with velocity increasing from 4.12 km/s to 4.57 km/s, where we interpret the Moho to be. Similar to other models, the Moho can either be gradational, where we don’t see a distinct layer with a huge jump in velocity, or it can be a sharp boundary that has been smoothened by the smoothness parameter of 0.70. The surface wave dispersion shows a misfit of 0.0724 for the Rayleigh waves and a 0.1442 for the Love waves. The higher misfit might be due to the modeled long periods not matching the data. The receiver functions show a misfit of 0.0406, where the Pp, Ps, and the reverberating phases’ arrival times fit but the amplitudes are lower for the final waveform. Something to note is that the resulting receiver function for a model using an “a priori” 66 Figure 32: Joint inversion result for group of stations WKZ05 and WKZ06. The box at left shows the initial a priori velocity model (gray) compared to the final velocity model (red). The green dashed lines represent the slope breaks that mark the boundary of the sections. From the final velocity model, a final receiver function (bottom right) and surface wave dispersion curves (top right) were calculated. The “initial” results are calculated from the a priori model to show what the results would be if the model was true. This final velocity model shows two upper sections of 5.5 km and 5 km, a 12.5 km thick middle section, a 10 km thick lower section, followed by a 7.5 km thick gradient at the base of the crust with an overall crustal thickness of 33 – 40 km. 67 weight of 0.3 decreases the misfit in the receiver function, but it creates more noise in the final velocity model. The second model, station WKZ07, is also inverted by a priori weight of 0.5, a weight of 50 – 50% to receiver functions and dispersion, and a smoothness of 0.70 (Figure 33). Two upper sections are seen: the top section starts at 1.79 km/s increasing to 2.40 km/s at 3 km depth and the lower section has a velocity of 2.88 km/s to 5.5 km depth. Then a middle section follows steadily increasing from 3.17 km/s to 3.88 km/s at 25.5 km depth. The lower section is more difficult to identify since the slope break is very subtle, but it is found from 25.5 to 35.5 km depth with velocity 3.95 km/s to 4.08 km/s. The base of the crust has a velocity gradient from 4.14 km/s at 36 km to 4.58 km/s at 45.5 km depth, which is where we interpret the Moho is. Since the receiver function stack from this station does not show the PpPs and PsSs+PsPs phases, the velocity model does not show a sharp Moho depth. The misfits for the Surface wave dispersion curves are 0.0915 for Rayleigh waves and 0.1437 for Love waves. The receiver function misfit is 0.0377, where the final receiver function has a lower Pp amplitude than the observed and the troughs between the Ps and reverberations and noise have a higher amplitude. After running multiple inversions with different parameters, we found that many things can affect the inversion result. Most of the groups produced reliable results, although they include error emerging from the scarcity of data used and the surface wave dispersion curve representing sedimentary basin velocities and not including shorter period waves. In general, shorter period surface waves can sample fine details in the upper Qualitative effect of data quality on results 68 Figure 33: Joint inversion result for station WKZ07. The box at left shows the initial a priori velocity model (gray) compared to the final velocity model (red). The green dashed lines represent the slope breaks that mark the boundary of the sections. From the final velocity model, a final receiver function (bottom right) and surface wave dispersion curves (top right) were calculated. The “initial” results are calculated from the a priori model to show what the results would be if the model was true. This final velocity model shows two upper sections of 3 km and 2.5 km, a 20 km thick middle section, a 10 km thick lower section, followed by a 10 km thick gradient at the base of the crust with an overall crustal thickness of 36 – 45 km. 69 crust but are not reliable in the surface wave dispersion dataset from Pasyanos (2005), while longer period waves sample regions larger than some structures we are trying to image. If we want a lower misfit of the dispersion curves, we need to make a model that produces a dispersion curve resembling the observed data, even though the model can be unrealistic. Nonetheless, the joint inversion provided us with very useful information about the crustal structure and shear wave velocity. The data used to process the receiver functions are very clean and quiet, but the process needs more data. Since we use a Gaussian filter of 1.0, receiver functions don’t resolve the shallower structure either, hence the reason for inverting with surface wave dispersion. The receiver functions contain their own noise and can add error for the processing of shear wave velocity models. The group of stations WKZ03 and WKZ04 were not interpreted since it provided unstable results and did not produce a coherent shear wave velocity model. This can be because the northwest part of Central Mangyshlak is where the crust is the most heterogenous according to Dimakov (1968) and the Mstislavskiy et al. (1980) cross section (see Figure 4). If we wanted to map the crust where stations WKZ03 and WKZ04 were located, we will need more data to produce receiver functions and get information on finer structure. The inversion for the group of WKZ05 and WKZ06 provided good results, but the misfit between the dispersion curves is a bit high. If we use an a priori weight of 0.3, the dispersion misfit decreases, but the velocity gradient at the base of the crust is broader. The trouble in this group is not with weights given to each dataset or smoothness, but with the receiver functions. 70 The stack of station WKZ07does not show well defined PpPs and PsSs+PsPs phases, which provide additional constraints on the Moho depth. Even though we had some errors in the datasets (i.e. missing phases, low amplitude phases for RF, or no short period surface waves), we were able to compute good results which can be compared to previous knowledge of this region. 71 We compare our results to previous studies for each area in our region in regard to crustal structure, seismic velocities, Moho depth and composition, and discuss factors that might have affected our results. We produced shear wave velocity structure models of the crust for each of our areas in the Mangystau Region. From previous studies we infer that these areas vary in structure and seismic velocity. We interpret the type of lithology at different depths of the crust based on their seismic velocities from the joint inversion of receiver functions with surface wave dispersion. The crustal structure was interpreted by separating different sections in the final velocity models where layer boundaries were determined to be at slope breaks, or inflection points. The velocities used in the processing of joint inversions resembled those of sedimentary and metamorphosed-sedimentary, “granitic”, and “basaltic” rocks. From this comparison of values, we created crustal columns that represent thickness and average velocity of each section within the crust and the rock type most likely to be present in these sections (Figure 34). The velocity structure results for North Ustyurt provide a reasonable model, almost resembling the a priori model perfectly. We compare these results with the seismic profile Aral-1 from Figure 6 (Figure 35) and the northern part of Dimakov (1968) Profile I from Figure 4 (Figure 36). From our results, we interpret that the crust in North Ustyurt consists of a 3 km sedimentary basin, on top of a 7 km thick Permian-Triassic Complex, a 10 km thick granitic mid-crust, and a 16 to 21 km thick basaltic lower crust, with an overall crustal thickness of 36 – 41 km. The sedimentary basin is about the same thickness as the North Ustyurt DISCUSSION OF RESULTS 72 Figure 34: Columns representing the crustal structure of the Mangystau Region from results. The crust in all 3 areas is divided into 5 sections: from top to bottom, a Jurassic to Quaternary Sedimentary basin, a Permian-Triassic metamorphosed sedimentary Complex, mid-crustal “granitic” upper basement, and “basaltic” lower basement. The types of rocks were determined from previously assumed lithologies in this region from Dimakov (1968). The black numbers inside the blocks are lower and upper Vs values, and the green numbers are the average layer Vs obtained from the joint inversion results. 73 Figure 35: Comparison of a crustal column from our results in North Ustyurt and a cross-section from Sheikh-Zade (1996). Both, column and cross-section, are in the same vertical scale. The column at left is the interpretation of our results and the cross- section to the right is from Sheikh-Zade (1996) showing a seismic study done in the Aral Sea. The black numbers in the right cross section are Vp measurements at depth, while the red numbers are Vs calculated with a Vp/Vs of 1.7. For legend of the cross-section see Figure 6. 74 Figure 36: Comparison of crustal column from our results of North Ustyurt and Dimakov (1968) Profile I. The left column is the interpretation from our results and is compared to the previous studies from Dimakov (1968) (right), showing differences in the structure of the crust, in the same vertical scale. The dashed line between the columns is used to help guide where the same boundary is relative to one another. We calculated a sedimentary basin and Permian-Triassic Complex very similar to Dimakov (1968). The difference in thickness can be due to the blurry image from the profile, but that ~0.5 km can’t be resolved with the smoothness we applied. Our granitic section is significantly thinner with a difference of about 6 km, but out basaltic section is almost the same, but we have an overall crustal thickness that is 5 – 9 km thinner, if we consider the whole gradient at the base of the crust. 75 Aral-1 profile and Profile I from Dimakov (1968). To the east of our stations, Sheikh-Zade (1996) did not find any internal structure in the basement that we can compare to, but we interpreted a boundary between the granitic mid-crust and basaltic lower crust at 20 km. To the west, Dimakov (1968) found a boundary between the granitic and basaltic units that is 5 km deeper than our interpreted boundary. Our calculated velocities for most of the crust are much slower than those from the Aral-1 profile. One reason might be that the rock type changes from our area to the Aral Sea; although Paleozoic sediments are mentioned to be present in the Aral Sea (Zheikh- Zade, 1966) these might not be the same that form the Permian – Triassic Complex. Our velocities, however, are very similar to those from Dimakov (1968), except for the basaltic lower basement, where we have a slower section. This might be caused by the velocity gradient at the base of the crust, where the Moho depth is not well constrained leading to miscalculating the depth and average velocity. One major difference is the crustal thickness, where we couldn’t resolve a Moho depth with the smoothing parameter we were using, while Dimakov (1968) reported a sharp 45 km deep Moho and Sheikh-Zade (1996) reported a gradational Moho at 35-37 km, although this was further east to the Aral Sea. Our model did not show a sharp velocity contrast at the Moho because of the smoothing parameter we used. When decreasing the smoothing parameter, many low-velocity layers appear in the middle of the crust, (Figure 14), and the Moho has a sharp velocity contrast. Previously mentioned studies show us that there are no such low-velocity layers at mid-crustal depths. We consider these layers to be artifacts, and in order to get rid of them we increased the smoothing parameter. The 76 Central Mangyshlak increased smoothness smears the Moho, thus we interpret the boundaries using slope changes. This is also true for Central Mangyshlak and South Mangyshlak. The receiver functions show different PpPs and PsSs+PsPs arrivals at different times (16 s and others at 17-18 s), which affect the amplitude of those phases in the stack. The phases are still present, but if they arrived at the same time in all the receiver functions, the stacks should be more coherent, and the amplitude might have been higher and might show the clear contrast from the Moho in the final velocity model. The reason for this arrival time discrepancy does not seem to be directional after comparing the arrival time of these phases of the same events for both stations WKZ09 and WKZ10, but more data is necessary to find evidence or the reason why these are shifted. Our results from Central Mangyshlak match the previous known information from Dimakov (1968) (Figure 37). The most representative model of the data was the grouping of WKZ01 and WKZ02. The stations were installed south of the uplifted segment in Central Mangyshlak, not directly over the thickest part, but they are between the two major faults that bound the Central Mangyshlak. Our model shows that the southern part of Central Mangyshlak consists of a 3 km thick upper sedimentary layer, over a 7.5 km Permian- Triassic Complex, a 12.5 km thick granitic mid crust, and 12.5 to 20 km basaltic lower crust, with an overall crustal thickness of 36 – 43 km. The Moho (Vs = 4.4 km/s) appears at a depth of 38 km, but it can be in the range of 36 to 43 km. Since we applied the smoothing parameter of 0.7, we can’t know exactly what the depth of the Moho is. The crustal thickness that we estimated is within the range of what Dimakov (1968) and Volvovsky et al. (1966) have assumed. 77 Figure 37: Comparison of crustal column from our results of Central Mangyshlak and Dimakov (1968) Profile I. The left column is the interpretation from our results and is compared to the previous studies from Dimakov (1968) (right), showing differences in the structure of the crust, in the same vertical scale. The dashed line between the columns is used to help guide where the same boundary is relative to one another. We calculated a sedimentary basin 1 km thicker, and Permian-Triassic Complex 2 km thicker to Dimakov (1968). This difference is not as significant since these boundaries were selected by slope breaks instead of sharp steps in the final velocity models. The basaltic sections are about the same thickness, but ours is 3 km deeper than Dimakov (1968). Our basaltic section, above the velocity gradient at the bottom of the crust, is thinner than Dimakov (1968), but if we assume that the bottom of the gradient is where the basaltic section ends, we would have a very similar thickness of this section. 78 Our average velocities are very similar to Dimakov (1968), except for the sedimentary basin. If we assumed the boundary of the sedimentary basin and the Permian- Triassic Complex to be shallower in the crust, we might be able to obtain an average velocity similar to Dimakov (1968), but then the Permian-Triassic velocity would be different. The same case is with the basaltic section, where we can get a higher velocity section, but we are not sure where this section ends. If we assume that the Moho where Vs=4.4 km/s, then the basaltic section might be similar to Dimakov (1968). Something interesting to point out is that stations WKZ03 and WKZ04 did not provide good results in this study. The receiver function stack from stations WKZ01 and WKZ02 and the stack of WKZ03and WKZ04 look very similar, as does the surface wave dispersion for both areas, but they produced very different shear wave velocity structures. The difference between the observed and final receiver functions and surface wave dispersion curves of the stack of stations WKZ03 and WKZ04 was high with all a priori models, if the models converged. The reason for this might be due to the crust under WKZ03 and WKZ04 being the most deformed in the whole region, causing uncertainties in the models. The shear wave velocity models for South Mangyshlak from the grouping of WKZ05 and WKZ06 and the one from WKZ07 converge to something reasonable and match somewhat with Dimakov (1968) Profile I. Based on the WKZ05 and WKZ06 results, South Mangyshlak consists of a 5.5km thick upper sedimentary basin, over a 5 km thick Permian- Triassic Complex, a 12.5 granitic section, and a 10 to 22.5 km thick basaltic section, with an overall crustal thickness of 33 – 40 km (Figure 38). The final velocity model of WK07 shows South Mangyshlak 79 Figure 38: Comparison of crustal columns from our results of South Mangyshlak and Dimakov (1968) Profile I. The left two columns are the interpretation from our results and are compared to the previous studies from Dimakov (1968) (right), showing differences in the structure of the crust, in the same vertical scale. The dashed line between the columns is used to help guide where the same boundary is relative to one another. From the group of WKZ05 and WKZ06, we calculated a sedimentary basin about 3 km thicker than Dimakov (1968), and a Permian-Triassic Complex similar in thickness, but it extends 4 km deeper. Our granitic section about 2 km thinner, but going down to the same depth, and a basaltic section with about the same thickness, if we consider the basaltic section going down to the bottom extent of the velocity contrast. Even with this assumption, the Moho from Dimakov (1968) is about 2 km thicker than our interpreted crustal thickness. The interpretation from WKZ07 shows a sedimentary basin with the same thickness as Dimakov (1968), a 1 km thinner Permian-Triassic Complex, a 4-5 km thicker granitic section, and a basaltic section of very similar depth. 80 a 3 km thick upper sedimentary basin, over a 2.5 km thick Permian-Triassic Complex, a 20 km thick granitic section, and a 10 to 20 km thick basaltic section, with an overall crustal thickness of 36 – 45 km (Figure 38). Some differences exist between our two models for this area. The sedimentary basin for both models share the same seismic velocity, but the velocity of the Permian-Triassic Complex is lower in WKZ07. Although the granitic sections from both of our columns share the same seismic velocity with Dimakov (1968), they vary in thickness. The granitic section from WKZ05-WKZ06 is very similar to the one from Dimakov (1968), but the one from WKZ07 is thicker than both of the former columns. Station WKZ07 is further southeast of both of Dimakov (1968) profiles, and stations WKZ05 and WKZ06. This difference in section thicknesses can represent different structure under station WKZ07. The basaltic section is comparable in size and velocity to the one from Dimakov (1968), but still has uncertainty in Moho depth due to the velocity gradient at the bottom of the crust. To better resolve the Moho transition, we need more receiver functions to be able to find a representative waveform for each station. In the case for all of all of our geomorphological areas, if we apply a low smoothness parameter value to the final velocity models, we generate low velocity layers in the mid crust, and if we increase it, we smear all the velocity layers and make the Moho look gradational (see Figure 14). Most of the individual receiver functions contain clear arrivals from the phases converted at the Moho, which implies a very distinct Moho. Once we stack the receiver functions to increase the signal to noise ratio, however, some of the phases interfere with each other and might decrease in amplitude making the signal seem like there is no distinct Moho. At least this is Moho discontinuity 81 the case for the grouping of stations WKZ05-WKZ06 and station WKZ07. On the other hand, a gradational boundary has already been imaged in the Aral Sea. Therefore, we conclude that with the amount of data used for this study, the Moho discontinuity depth and “sharpness” in the Mangystau Region cannot be imaged precisely. 82 We were able to determine the crustal structure of the Mangystau Region in southwestern Kazakhstan by computing receiver functions and inverting them jointly with surface wave dispersion curves, using the method of Julià et al. (2000), and comparing the results with previous studies from Sheikh-Zade (1996) and Dimakov (1968), among others. The seismic data available for the calculation of receiver functions was limited but it was still useful since it was exceptionally quiet. Surface wave dispersion curves were obtained from Pasyanos (2005) and although the values seem reliable and help converge to reasonable velocity modes, there weren’t any stations in the region that were used to calculate such values. The surface wave dispersion curves we obtained were interpolated from stations surrounding the region, but they were still useful. In order to overcome this issue of data availability and dispersion curves, we need to collect more data in the region and calculate surface wave dispersion curves with such data. The smoothness in the final velocity models had to be increased due to noise in the waveforms after stacking but can be decreased if more interpretative receiver functions are computed to extract more detail. Also, with more data, we will be able to compute receiver functions and shear wave velocity models for each individual station instead of having to group some together. Although the crustal structure was previously imaged during the 1960’s, we applied newer methods to see if these are useful in the region and see if there are different interpretations. Our results were very similar to previously calculated ones, and we sampled further north in North Ustyurt than Dimakov (1968). We were able to determine layer depths and seismic velocities within the crust to some degree. Surprisingly, the CONCLUSIONS 83 shallower layers were more consistent with previous studies than the deeper ones, even though we did not have the 5-15 second period surface waves that sample such depths. All models matched the Dimakov (1968) profile to a certain extent but were not completely the same. Since we were marking the boundaries where there were slope breaks, it was difficult to define where these were the slope breaks were not so sudden, allowing us to fit the velocities and depths that we thought were appropriate. The one thing that we were not able to interpret confidently was the Moho depth. This was due to the smoothness parameter applied to decrease the lower velocity layers created in the mid-crust, and because of the receiver functions stacks. We cannot say whether the Moho has a gradational nature of it is a sharp boundary at the base of the crust. We showed that joint inversions of receiver functions and surface wave dispersion can help us find information of the crustal structure of a region. Fortunately, we had clean, high-quality data that we could process to do this study, but we still did not have enough to interpret a Moho boundary or more distinct layers. With more data, and with the same or better quality, we would be able to conclude whether our models are more accurate, or if they had errors and Dimakov (1968) has the better results. Since the seismic data collected was very quiet and such good signal, a new deployment in the area is planned for 2020. This deployment will consist of installing 16 broadband seismic stations for one year, which will allow us to investigate the crustal structure further. Additionally, we plan on studying if there is any anisotropy in the upper mantle via shear wave splitting to infer the direction of maximum finite strength. This 84 might give us an indication of what has happened in the past tectonically, since there is no active tectonic activity occurring in this region. 85 BIBLIOGRAPHY 86 BIBLIOGRAPHY Abetov A., G. Zhylkybayeva, and T. Zhylkybayev (2015). Deep structure, seismic models, rheology and geodynamics of consolidated crust. Canadian Scientific Journal 1, 65 – 77. Ammon, C. J., G. E. Randall, and G. Zandt (1990). On the nonuniqueness of receiver function inversions. Journal of Geophysical Research 95, 15303 – 15318. doi:10.1029/JB095iB10p15303 Babazhanov T. L., V. A. Rzaeva, and E. R. Sheikh-Zade (1989). Glubinnoe stroenie zemloi kori pod Aral’skim morem. Sovetskaya Geologiya (1) 80 – 86. (in Russian). Babadjanov T. L., V. A. Rzaeva, M. M. Rzaeb, E. R. Shiekh-Zade (1991). Metod glubinnikh otrazhennikh voln pri nzuchenni zemnoi kori Vostochnogo Ustyurta. Sovetskaya Geologiya (3) 56 – 63. (in Russian). Bazhenov, M. L., V. S. Burtman, and A. V. Dvorova (1999). Permian paleomagnetism of the Tien Shan fold belt, Central Asia: post-collisional rotations and deformation. Tectonophysics 312, 303 – 329. doi:10.1016/S0040-1951(99)00181-X Bird, P. (2003). An updated digital model of plate boundaries. Geochemistry, Geophysics, Geosystems 4. doi:10.1029/2001GC000252 Chai, C., C. J. Ammon, M. Maceira, and R. B. Herrmann (2015). Inverting interpolated receiver functions with surface wave dispersion and gravity: Application to the western U.S. and adjacent Canada and Mexico. Geophysical Research Letters 42, 4359 – 4366. doi:10.1002/2015GL063733 Dimakov, A. I. (1968). Stroenie zemnoi kory v raione Mangyshlaka. Sovetskaya Geologiya (11) 122 – 127 (In Russian). Dziewonski, A. M., and D. L. Anderson (1981). Preliminary reference Earth Model. Physics of the Earth and Planetary Interiors 25, 297 – 356. doi:10.1016/0031-9201(81)90046- 7 Gaetani, M., M. Balini, V. J. Vuks, V. A. Gavrilova, E. Garzanti, A. Nicora, E. Erba, E. Cariou, F. Cecca, S. I. Premoli, M. R. Pet- rizzo, S. Cirilli, and B. P. Raffaella (1998). The Mesozoic of the Mangyshlak (West Kazakhstan). Mémoires du Muséum national d'histoire naturelle 179, 35-74. Golonka, J. (2004). Plate tectonic evolution of the southern margin of Eurasia in the Mesozoic and Cenozoic. Tectonophysics 381, 235 – 273. doi:10.1016/j.tecto.2002.06.004 87 Julià, J., C. J. Ammon, R. B. Herrmann, and A. M. Correig (2000). Joint inversion of receiver function and surface wave dispersion observations. Geophysical Journal International 143, 99 -112. doi:10.1046/j.1365-246x.2000.00217.x Khain, V. E. (1994). Geology of Northern Eurasia (Ex-USSR). Gebruder Berntraeger, Berlin, 404 p. Kikuchi, M., and H. Kanamori (1982). Inversion of complex body waves. Bulletin of Seismological Society of America 72, 491 – 506. Langston, C. A. (1979). Structure under Mount Rainier, Washington, inferred from Teleseismic Body Waves. Journal of Geophysical Research 84, 4749 – 4762. doi:10.1029/JB084iB09p04749 Lemarie, M. M., M. Westphal, E. L. Gurevitch, K. Nazarov, H. Feinberg, and J. P. Pozzi (1997). How far between Iran and Eurasia was the Turan plate during Triassic–Jurassic times? Geologie en Mijnbouw 76, 73 – 82. doi:10.1023/A:1003129028277 Ligorría, J. P. and C. J. Ammon (1999). Iterative Deconvolution and Receiver-Function Estimation. Bulletin of the Seismological Society of America 89, 1395 – 1400. Mackey, K., I. Sokolova, D. Burk, A. Abishev, A. Belyashov, and R. Jih (2017). A Seismic Noise Survey of Western Kazakhstan. CTBTO: Science and Technology Conference 2017. T3.1-P5. Maceira, M., and C. J. Ammon (2009). Joint inversion of surface wave velocity and gravity observations and its application to central Asian basins shear velocity structure. Journal of Geophysical Research 114. doi:10.1029/2007JB005157 Mstislavskiy, M. M., A. M. Mezentsev, A. M., and L. N. Olofinskiy (1980). Tectonic regionalisation of the pre-jurassic basement in the Mangyshlak uplift zone. Geotectonics 14, 121 – 130. Natal’in, B. A., and A. M. Celâl Şengör (2005). Late Palaeozoic to Triassic evolution of the Turan and Scythian platforms: The pre-history of the Palaeo-Tethyan closure. Tectonophysics 404, 175 – 202. doi:10.1016/j.tecto.2005.04.011 Pasyanos, M. E., W. R. Walter, and S. E. Hazler (2001). A surface wave dispersion study of the Middle East and North Africa for monitoring the Comprehensive Nuclear-Test- Ban Treaty. Pure and Applied Geophysics 158, 1445 – 1474. doi:10.1007/PL00001229 Pasyanos, M. E. (2005). A variable resolution surface wave dispersion study of Eurasia, North Africa, and surrounding regions. Journal of Geophysical Research 110. doi:10.1029/2005JB003749 88 Popkov, V. I. (1986). Tectonics of the Pre-Jurassic Sedimentary Complex in the Western Turanian Plate. Geotectonics 20, 328 – 335. Popkov, V. I. (1995). Tectonics of the west of the Turan platform. Petroleum Geology 29, 71 – 102. Randall, G. E. (1994). Efficient Calculation of complete differential seismograms for laterally homogeneous Earth models. Geophysical Journal International 118, 245 – 254. doi:10.1111/j.1365-246X.1994.tb04687.x Ryaboy, V. Z. (1967). Structure of Earth’s crust and upper mantle along deep seismic profile (Kopet-Dag to Aral Sea). International Geology Review 9, 296 – 299. doi:10.1080/00206816709474467 Ryaboy V. Z. (1989). Upper mantle structure studies by explosion seismology in the USSR. Delphic Assoiciates, Inc., Falls Church Virginia, USA, 141 p. Sheikh-Zade, E. R. (1996). Results of seismic reflection profiling in the Turarnian Platform. Tectonophysics 264, 123 – 135. doi:10.1016/S0040-1951(96)00122-9 Shen, W., M. H. Ritzwoller, V. Schulte-Pelkum, and F. Lin (2013). Joint inversion of surface wave dispersion and receiver functions: a Bayesian Monte-Carlo approach. Geophysical Journal International 192, 807 – 836. doi:10.1093/gji/ggs050 Thomas, J. C., P. R. Cobbold, V. S. Shein, and S. L. Douaran (1999). Sedimentary record of late Paleozoic to Recent tectonism in central Asia - analysis of subsurface data from the Turan and south Kazak domains. Tectonophysics 313, 243 – 263. doi:10.1016/S0040-1951(99)00208-5 USGS NEIC (accessed September 2016). Search Earthquake Catalog. Retrieved from https://earthquake.usgs.gov/earthquakes/search/ Volvovsky, I. S., R. G. Garetzkiy, A. E. Shlezinger, and V. N. Shraybman (1966). Tektonika Turanskoy Plity. Akademiya Nauka SSSR, Geologicheskiy Institut Trudy 165, 288 p, “Nauka” Publishers, Moscow, USSR (In Russian). Tectonics of The Turan Plate; Technical Translation by Synenko G., and D. Shalikashvili (1969). Linguistic Section, Western Area Branch, Chart Research Division, Aeronautical Chart and Information Center, St. Louis, Missouri, USA, 392 p. Wessel, P., W. H. F. Smith, R. Scharroo, J. F. Luis, and F. Wobbe (2013). Generic Mapping Tools: Improved version released. EOS Transactions, American Geophysical Union 94, 409-410. doi:10.1002/2013EO450001 Yuferov, Y. K., V. Y. Aronson, and A. A. Rabinovich (1973). Distribution of oil and gas pools in the Zhetybay-Uzen zone of oil-gas accumulation. Petroleum Geology 11, 197 – 201. 89 Zonenshain, L. P., M. I. Kuzmin, and L. M. Natapov (1990). Geology of the USSR: a plate tectonic synthesis. AGU Geodynamic Series 21, American Geophysical Union, Washington DC, 242 p. doi: 10.1029/GD021 Zuza, A. V., and A. Yin (2017). Balkatach hypothesis: A new model for the evolution of the Pacific, Tethyan, and Paleo-Asian oceanic domains. Geosphere 13, 1664 – 1712. doi:10.1130/GES01463.1 90