ALBERT

Description: Trade-offs between velocity and anisotropy heterogeneity complicate the interpretation of differential traveltime data and have the potential to bias isotropic tomographic models. By constructing a simple parametrisation to describe an elastic tensor with hexagonal symmetry, we find analytic solutions to the Christoffel equations in terms of fast and slow horizontal velocities that allow us to simultaneously invert differential traveltime data and splitting data from teleseismic S arrivals to recover 3-D velocity and anisotropy structure. This technique provides a constraint on the depth-extent of shallow anisotropy, otherwise absent from interpretations based on SKS splitting alone. This approach is well suited to the young Woodlark Rift, where previous studies have found strong velocity variation and substantial SKS splitting in a continental rift with relatively simple geometry. This study images a low-velocity rift axis with ≤4 per cent spreading-parallel anisotropy at 50–100 km depth that separates regions of pre-existing lithospheric fabric, indicating the synchronous development of extensional crystallographic preferred orientation and lithospheric thinning. A high-velocity slab fragment north of the rift axis is associated with strike-parallel anisotropic fast axes, similar to that seen in the shallow mantle of some subduction zones. In addition to the insights provided by the anisotropy structure, the improvement in fit to the differential traveltime data demonstrates the merit to a joint inversion that accounts for anisotropy.

Description: Short-period small magnitude seismograms mainly comprise scattered waves in the form of coda waves (the tail part of the seismogram, starting after S waves and ending when the noise prevails), spanning more than 70 per cent of the whole seismogram duration. Corresponding coda envelopes provide important information about the earth inhomogeneity, which can be stochastically modeled in terms of distribution of scatterers in a random medium. In suitable experimental conditions (i.e. high earth heterogeneity), either the two parameters describing heterogeneity (scattering coefficient), intrinsic energy dissipation (coefficient of intrinsic attenuation) or a combination of them (extinction length and seismic albedo) can be used to image Earth structures. Once a set of such parameter couples has been measured in a given area and for a number of sources and receivers, imaging their space distribution with standard methods is straightforward. However, as for finite-frequency and full-waveform tomography, the essential problem for a correct imaging is the determination of the weighting function describing the spatial sensitivity of observable data to scattering and absorption anomalies. Due to the nature of coda waves, the measured parameter couple can be seen as a weighted space average of the real parameters characterizing the rock volumes illuminated by the scattered waves. This paper uses the Monte Carlo numerical solution of the Energy Transport Equation to find approximate but realistic 2-D space-weighting functions for coda waves. Separate images for scattering and absorption based on these sensitivity functions are then compared with those obtained with commonly used sensitivity functions in an application to data from an active seismic experiment carried out at Deception Island (Antarctica). Results show that these novel functions are based on a reliable and physically grounded method to image magnitude and shape of scattering and absorption anomalies. Their extension to 3-D holds promise to improve our ability to model volcanic structures using coda waves.

Description: We investigate the relationship between seismic moment M 0 and source duration t w of microearthquakes by using high-quality seismic data recorded with a vertical borehole array installed in central Taiwan. We apply a waveform cross-correlation method to the three-component records and identify several event clusters with high waveform similarity, with event magnitudes ranging from 0.3 to 2.0. Three clusters—Clusters A, B and C—contain 11, 8 and 6 events with similar waveforms, respectively. To determine how M 0 scales with t w , we remove path effects by using a path-averaged Q . The results indicate a nearly constant t w for events within each cluster, regardless of M 0 , with mean values of t w being 0.058, 0.056 and 0.034 s for Clusters A, B and C, respectively. Constant t w , independent of M 0 , violates the commonly used scaling relation ${t_w} \propto M_0^{1/3}$ . This constant duration may arise either because all events in a cluster are hosted on the same isolated seismogenic patch, or because the events are driven by external factors of constant duration, such as fluid injections into the fault zone. It may also be related to the earthquake nucleation size.

Description: Earthquake absolute location errors which can be encountered in an underground reservoir are investigated. In such an exploitation context, earthquake hypocentre errors can have an impact on the field development and economic consequences. The approach using the state-of-the-art techniques covers both the location uncertainty and the location inaccuracy—or bias—problematics. It consists, first, in creating a 3-D synthetic seismic cloud of events in the reservoir and calculating the seismic traveltimes to a monitoring network assuming certain propagation conditions. In a second phase, the earthquakes are relocated with assumptions different from the initial conditions. Finally, the initial and relocated hypocentres are compared. As a result, location errors driven by the seismic onset time picking uncertainties and inaccuracies are quantified in 3-D. Effects induced by erroneous assumptions associated with the velocity model are also modelled. In particular, 1-D velocity model uncertainties, a local 3-D perturbation of the velocity and a 3-D geostructural model are considered. The present approach is applied to the site of Rittershoffen (Alsace, France), which is one of the deep geothermal fields existing in the Upper Rhine Graben. This example allows setting realistic scenarios based on the knowledge of the site. In that case, the zone of interest, monitored by an existing seismic network, ranges between 1 and 5 km depth in a radius of 2 km around a geothermal well. Well log data provided a reference 1-D velocity model used for the synthetic earthquake relocation. The 3-D analysis highlights the role played by the seismic network coverage and the velocity model in the amplitude and orientation of the location uncertainties and inaccuracies at subsurface levels. The location errors are neither isotropic nor aleatoric in the zone of interest. This suggests that although location inaccuracies may be smaller than location uncertainties, both quantities can have a cumulative effect. Besides, small velocity uncertainties applied to the whole 1-D profile can lead to large increase of the location uncertainties. However, local variations of the velocity field around the well may have negligible effects that would make such a feature undetectable with an absolute location method. Although the reference 1-D velocity model was built from well log data, the results show that it is not a good representative of a more realistic 3-D model including a fault and its associated block shift. The amplitude and distribution of the induced location inaccuracies are such that the positioning and the orientation of features delineated by seismicity are distorted and may be difficult to correctly interpret.

Description: Understanding telomere length maintenance mechanisms is central in cancer biology as their dysregulation is one of the hallmarks for immortalization of cancer cells. Important for this well-balanced control is the transcriptional regulation of the telomerase genes. We integrated Mixed Integer Linear Programming models into a comparative machine learning based approach to identify regulatory interactions that best explain the discrepancy of telomerase transcript levels in yeast mutants with deleted regulators showing aberrant telomere length, when compared to mutants with normal telomere length. We uncover novel regulators of telomerase expression, several of which affect histone levels or modifications. In particular, our results point to the transcription factors Sum1, Hst1 and Srb2 as being important for the regulation of EST1 transcription, and we validated the effect of Sum1 experimentally. We compiled our machine learning method leading to a user friendly package for R which can straightforwardly be applied to similar problems integrating gene regulator binding information and expression profiles of samples of e.g. different phenotypes, diseases or treatments.

Description: The ability to integrate ‘omics’ (i.e. transcriptomics and proteomics) is becoming increasingly important to the understanding of regulatory mechanisms. There are currently no tools available to identify differentially expressed genes (DEGs) across different ‘omics’ data types or multi-dimensional data including time courses. We present fCI (f-divergence Cut-out Index), a model capable of simultaneously identifying DEGs from continuous and discrete transcriptomic, proteomic and integrated proteogenomic data. We show that fCI can be used across multiple diverse sets of data and can unambiguously find genes that show functional modulation, developmental changes or misregulation. Applying fCI to several proteogenomics datasets, we identified a number of important genes that showed distinctive regulation patterns. The package fCI is available at R Bioconductor and http://software.steenlab.org/fCI/ .

Description: Next generation sequencing of cellular RNA is making it possible to characterize genes and alternative splicing in unprecedented detail. However, designing bioinformatics tools to accurately capture splicing variation has proven difficult. Current programs can find major isoforms of a gene but miss lower abundance variants, or are sensitive but imprecise. CLASS2 is a novel open source tool for accurate genome-guided transcriptome assembly from RNA-seq reads based on the model of splice graph. An extension of our program CLASS, CLASS2 jointly optimizes read patterns and the number of supporting reads to score and prioritize transcripts, implemented in a novel, scalable and efficient dynamic programming algorithm. When compared against reference programs, CLASS2 had the best overall accuracy and could detect up to twice as many splicing events with precision similar to the best reference program. Notably, it was the only tool to produce consistently reliable transcript models for a wide range of applications and sequencing strategies, including ribosomal RNA-depleted samples. Lightweight and multi-threaded, CLASS2 requires 〈3GB RAM and can analyze a 350 million read set within hours, and can be widely applied to transcriptomics studies ranging from clinical RNA sequencing, to alternative splicing analyses, and to the annotation of new genomes.

Description: Compared with a single GPS system, GPS/GLONASS observations can improve the satellite visibility, optimize the spatial geometry and improve the precise positioning performance. Although the advantage over GPS-only methods in terms of positioning is clear, the potential contributions of GPS/GLONASS to co-seismic displacement determination and the subsequent seismic source inversion still require extensive study and validation. In this paper, we first extended a temporal point positioning model from GPS-only to GPS/GLONASS observations. Using this new model, the performance of the GPS/GLONASS method for obtaining co-seismic displacements was then validated via eight outdoor experiments on a shaking table. Our result reveals that the GPS/GLONASS method provides more accurate and robust co-seismic displacements than the GPS-only observations in a non-optimal observation environment. Furthermore, as a case study, observation data recorded during the September 2015 M w 8.3 Illapel earthquake in Chile were re-processed. At some stations, obvious biases were found between the co-seismic displacements derived from GPS-only and GPS/GLONASS observations. The subsequent slip distribution inversion on a curved fault confirms that the differences in the co-seismic displacements causes differences in the inversion results and that the slip distributions of the Illapel earthquake inferred from the GPS/GLONASS observations tend to be shallower and larger.

Description: Incorporation of attenuation into the normal mode sum representations of seismic signals is commonly effected by applying perturbation theory. This is fine for weak attenuation, but problematic for stronger attenuation. In this work, modes of the anelastic medium are represented as complex superpositions of elastic eigenfunctions. For the P – SV system, a generalized eigenvalue equation for the complex eigenwavenumbers and complex coefficients used to construct the anelastic eigenfunctions is derived. The generalized eigenvalue problem for the P – SV problem is exactly linear in the eigenwavenumber at the expense of doubling the dimension. The SH problem is exactly linear in the square of the eigenwavenumber. This is in contrast to a similar standing wave problem for the earth free oscillations. Attenuation is commonly incorporated into synthetic seismogram calculations by introduction of complex frequency-dependent elastic moduli. The moduli depend nonlinearly on the frequency. The independent variable in the standing wave free oscillation problem is the frequency, which makes the eigenvalue problem nonlinear. The choice of the wavenumber as the independent variable for the travelling wave problem leads to a linear problem. The Earth model may be transversely isotropic. Compressional waves and both polarizations of shear waves ( SV, SH ) are treated.

Description: We present a new technique for deriving detailed information on seismic velocities of the subsurface material from continuous ambient noise recorded by spatially dense seismic arrays. This method uses iterative double beamforming between various subarrays to extract surface wave contributions from the ambient-noise data in complex environments with unfavourable noise-source distributions. The iterative double beamforming extraction makes it possible to retrieve large amounts of Rayleigh wave traveltime information in a wide frequency band. The method is applied to data recorded by a highly dense Nodal array with 1108 vertical geophones, centred on the damage zone of the Clark branch of the San Jacinto Fault Zone south of Anza, California. The array covers a region of ~650 x 700 m 2 , with instrument spacing of 10–30 m, and continuous recording at 500 samples s –1 over 30 d in 2014. Using this iterative double beamforming on subarrays of 25 sensors and cross-correlations between all of the station pairs, we separate surface waves from body waves that are abundant in the raw cross-correlation data. Focusing solely on surface waves, maps of traveltimes are obtained at different frequencies with unprecedented accuracy at each point of a 15-m-spacing grid. Group velocity inversions at 2–4 Hz reveal depth and lateral variations in the structural properties within and around the San Jacinto Fault Zone in the study area. This method can be used over wider frequency ranges and can be combined with other imaging techniques, such as eikonal tomography, to provide unprecedented detailed structural images of the subsurface material.

Description: Knowledge of the quality factor of near-surface materials is of fundamental interest in various applications. Attenuation can be very strong close to the surface and thus needs to be properly assessed. In recent years, several researchers have studied the retrieval of attenuation coefficients from the cross correlation of ambient seismic noise. Yet, the determination of exact amplitude information from noise-correlation functions is, in contrast to the extraction of traveltimes, not trivial. Most of the studies estimated attenuation coefficients on the regional scale and within the microseism band. In this paper, we investigate the possibility to derive attenuation coefficients from seismic noise at much shallower depths and higher frequencies (〉1 Hz). The Euroseistest area in northern Greece offers ideal conditions to study quality factor retrieval from ambient noise for different rock types. Correlations are computed between the stations of a small scale array experiment (station spacings 〈2 km) that was carried out in the Euroseistest area in 2011. We employ the correlation of the coda of the correlation (C 3 ) method instead of simple cross correlations to mitigate the effect of uneven noise source distributions on the correlation amplitude. Transient removal and temporal flattening are applied instead of 1-bit normalization in order to retain relative amplitudes. The C 3 method leads to improved correlation results (higher signal-to-noise ratio and improved time symmetry) compared to simple cross correlations. The C 3 functions are rotated from the ZNE to the ZRT system and we focus on Love wave arrivals on the transverse component and on Love wave quality factors Q L . The analysis is performed for selected stations being either situated on soft soil or on weathered rock. Phase slowness is extracted using a slant-stack method. Attenuation parameters are inferred by inspecting the relative amplitude decay of Love waves with increasing interstation distance. We observe that the attenuation coefficient and Q L can be reliably extracted for stations situated on soft soil whereas the derivation of attenuation parameters is more problematic for stations that are located on weathered rock. The results are in acceptable conformance with theoretical Love wave attenuation curves that were computed using 1-D shear wave velocity and quality factor profiles from the Euroseistest area.

Description: Large earthquakes from the intermediate-depth Vrancea seismic zone are known to produce in Bucharest ground motion characterized by predominant long periods. This phenomenon has been interpreted as the combined effect of both seismic source properties and site response of the large sedimentary basin. The thickness of the unconsolidated Quaternary deposits beneath the city is more than 200 m, the total depth of sediments is more than 1000 m. Complex basin geometry and the low seismic wave velocities of the sediments are primarily responsible for the large amplification and long duration experienced during earthquakes. For a better understanding of the geological structure under Bucharest, a number of investigations using non-invasive methods have been carried out. With the goal to analyse and extract the polarization and dispersion characteristics of the surface waves, ambient vibrations and low-magnitude earthquakes have been investigated using single station and array techniques. Love and Rayleigh dispersion curves (including higher modes), Rayleigh waves ellipticity and SH -wave fundamental frequency of resonance ( f 0 SH ) have been inverted simultaneously to estimate the shear wave velocity structure under Bucharest down to a depth of about 8 km. Information from existing borehole logs was used as prior to reduce the non-uniqueness of the inversion and to constrain the shallow part of the velocity model (〈300 m). In this study, we use data from a 35-km diameter array (the URS experiment) installed by the National Institute for Earth Physics and by the Karlsruhe Institute of Technology during 10 months in the period 2003–2004. The array consisted of 32 three-component seismological stations, deployed in the urban area of Bucharest and adjacent zones. The large size of the array and the broad-band nature of the available sensors gave us the possibility to characterize the surface wave dispersion at very low frequencies (0.05–1 Hz) using frequency–wavenumber techniques. This is essential to explore and resolve the deeper portions of the basin. The horizontal to vertical spectral ratio (H/V) curves provide important additional information about the structure and are here characterized by two major peaks. The first is attributed to the fundamental frequency of the basin, while the second can be interpreted as a mixture of the second higher mode of Rayleigh waves and other types of waves such as SH waves. This hypothesis has been verified by comparing the H/V curves with the SH -wave transfer function from the retrieved velocity structure. We could also approximate the SH transfer function with H/V ratios of earthquake recordings, providing additional verification of the robustness of the proposed velocity model. The Cretaceous bedrock depth was then inverted at each URS station from the fundamental frequency of resonance and using this model. A 3-D geophysical model for Bucharest has been constructed based on the integration of the inverted velocity profiles and the available geological information using a geographic information system.

Description: Some modern seismicity in the magnitude range of M 4 and M 6 in California and eastern North America preferentially occurs at the edges of past large ruptures. Once a large earthquake rupture has occurred, stress is concentrated at the edges of the rupture, and apparently this stress concentration can trigger earthquakes at or near the rupture edges many decades or even longer after a main shock. Furthermore, the modern M ≥ 4 earthquakes in the vicinity of a past main shock usually have the same focal mechanism as the earlier main shock. There are a number of examples of this in California and Nevada, where there is a statistically significant correlation of the locations of M ≥ 4 earthquakes and the edges of 19th and 20th century fault ruptures in M w ≥ 6.5 earthquakes. In contrast, the M ≥ 4 earthquakes near the epicentres of future ruptures in California are randomly scattered around the fault with no concentration near the ends of the future fault rupture. The concentration of earthquakes near the ends of earlier large ruptures in California becomes progressively less pronounced as the smallest magnitude in the data set is reduced from M 4.0 to M 3.0. These observations also appear to be true for intraplate regions where aftershock sequences can last millennia. The identification of modern rupture-edge M ≥ 4 aftershocks can be used to help discover where and when past strong earthquakes took place, even if there is no historical record of the main shock. This is of great importance for seismic hazard studies.

Description: Allele-specific copy number analysis (ASCN) from next generation sequencing (NGS) data can greatly extend the utility of NGS beyond the identification of mutations to precisely annotate the genome for the detection of homozygous/heterozygous deletions, copy-neutral loss-of-heterozygosity (LOH), allele-specific gains/amplifications. In addition, as targeted gene panels are increasingly used in clinical sequencing studies for the detection of ‘actionable’ mutations and copy number alterations to guide treatment decisions, accurate, tumor purity-, ploidy- and clonal heterogeneity-adjusted integer copy number calls are greatly needed to more reliably interpret NGS-based cancer gene copy number data in the context of clinical sequencing. We developed FACETS, an ASCN tool and open-source software with a broad application to whole genome, whole-exome, as well as targeted panel sequencing platforms. It is a fully integrated stand-alone pipeline that includes sequencing BAM file post-processing, joint segmentation of total- and allele-specific read counts, and integer copy number calls corrected for tumor purity, ploidy and clonal heterogeneity, with comprehensive output and integrated visualization. We demonstrate the application of FACETS using The Cancer Genome Atlas (TCGA) whole-exome sequencing of lung adenocarcinoma samples. We also demonstrate its application to a clinical sequencing platform based on a targeted gene panel.

Description: Surface wave methods provide a cost effective means of developing shear wave velocity ( Vs ) profiles for applications such as dynamic site characterization and seismic site response analyses. However, the inverse problem involved in obtaining a realistic layered earth model from surface wave dispersion data is inherently ill-posed, non-linear and mix-determined, without a unique solution. When available, a priori information such as geotechnical boreholes or geologic well logs should be used to aid in constraining site-specific inversion parameters. Unfortunately, a priori information is often unavailable, particularly at significant depths, and a ‘blind analysis’ must be performed. In these situations, the analyst must decide on an appropriate number of layers and ranges for their corresponding inversion parameters (i.e. trial number of layers and ranges in their respective thicknesses, shear wave velocities, compression wave velocities and mass densities). Selection of these parameters has been shown to significantly impact the results of an inversion. This paper presents a method for conducting multiple inversions utilizing systematically varied inversion layering parametrizations in order to identify and encompass the most reasonable layered earth models for a site. Each parametrization is defined by a unique layering ratio, which represents a multiplier that systemically increases the potential thickness of each layer in the inversion parametrization based on the potential thickness of the layer directly above it. The layering ratio method is demonstrated at two sites associated with the InterPacific Project, wherein it is shown to significantly aid in selecting reasonable Vs profiles that are close representations of the subsurface. While the goal of the layering ratio inversion methodology is not necessarily to find the ‘optimal’ or ‘best’ Vs profile for a site, it may be successful at doing so for certain sites/datasets. However, the primary reason for using the layering ratio method is to find Vs profiles that realistically represent the uncertainty in Vs resulting from surface wave inversion, and to avoid selection of Vs profiles that are unrealistic and adversely influenced by the choice of inversion parametrization.

Description: We present SWAN, a statistical framework for robust detection of genomic structural variants in next-generation sequencing data and an analysis of mid-range size insertion and deletions (〈10 Kb) for whole genome analysis and DNA mixtures. To identify these mid-range size events, SWAN collectively uses information from read-pair, read-depth and one end mapped reads through statistical likelihoods based on Poisson field models. SWAN also uses soft-clip/split read remapping to supplement the likelihood analysis and determine variant boundaries. The accuracy of SWAN is demonstrated by in silico spike-ins and by identification of known variants in the NA12878 genome. We used SWAN to identify a series of novel set of mid-range insertion/deletion detection that were confirmed by targeted deep re-sequencing. An R package implementation of SWAN is open source and freely available.

Description: The appraisal of tomographic models, of fundamental importance towards better understanding the Earth's interior, consists in analysing their resolution and covariance. The discrete theory of Backus–Gilbert, solving all at once the linear problems of model estimation and appraisal, aims at evaluating weighted averages of the true model parameters. Contrary to damped least-squares techniques, one key advantage of Backus–Gilbert inversion is that no subjective regularization is needed to remove the non-uniqueness of the model solution. Indeed, it is often possible to identify unique linear combinations of the parameters even when the parameters themselves are not uniquely defined. In other words, the non-uniqueness can be broken by averaging rather than regularizing. Over the past few decades, many authors have considered that, in addition to a high computational cost, it could be a clumsy affair in the presence of data errors to practically implement the Backus–Gilbert approach to large-scale tomographic applications. In this study, we introduce and adapt to seismic tomography the Subtractive Optimally Localized Averages (SOLA) method, an alternative Backus–Gilbert formulation which retains all its advantages, but is more computationally efficient and versatile in the explicit construction of averaging kernels. As a leitmotiv, we focus on global-scale S -wave tomography and show that the SOLA method can successfully be applied to large-scale, linear and discrete tomographic problems.

Description: High-throughput screening (HTS) is an indispensable tool for drug (target) discovery that currently lacks user-friendly software tools for the robust identification of putative hits from HTS experiments and for the interpretation of these findings in the context of systems biology. We developed HiTSeekR as a one-stop solution for chemical compound screens, siRNA knock-down and CRISPR/Cas9 knock-out screens, as well as microRNA inhibitor and -mimics screens. We chose three use cases that demonstrate the potential of HiTSeekR to fully exploit HTS screening data in quite heterogeneous contexts to generate novel hypotheses for follow-up experiments: (i) a genome-wide RNAi screen to uncover modulators of TNFα, (ii) a combined siRNA and miRNA mimics screen on vorinostat resistance and (iii) a small compound screen on KRAS synthetic lethality. HiTSeekR is publicly available at http://hitseekr.compbio.sdu.dk . It is the first approach to close the gap between raw data processing, network enrichment and wet lab target generation for various HTS screen types.

Description: The envelope broadening and the peak delay of the S -wavelet of a small earthquake with increasing travel distance are results of scattering by random velocity inhomogeneities in the earth medium. As a simple mathematical model, Sato proposed a new stochastic synthesis of the scalar wavelet envelope in 3-D von Kármán type random media when the centre wavenumber of the wavelet is in the power-law spectral range of the random velocity fluctuation. The essential idea is to split the random medium spectrum into two components using the centre wavenumber as a reference: the long-scale (low-wavenumber spectral) component produces the peak delay and the envelope broadening by multiple scattering around the forward direction; the short-scale (high-wavenumber spectral) component attenuates wave amplitude by wide angle scattering. The former is calculated by the Markov approximation based on the parabolic approximation and the latter is calculated by the Born approximation. Here, we extend the theory for the envelope synthesis of a wavelet in 2-D random media, which makes it easy to compare with finite difference (FD) simulation results. The synthetic wavelet envelope is analytically written by using the random medium parameters in the angular frequency domain. For the case that the power spectral density function of the random velocity fluctuation has a steep roll-off at large wavenumbers, the envelope broadening is small and frequency independent, and scattering attenuation is weak. For the case of a small roll-off, however, the envelope broadening is large and increases with frequency, and the scattering attenuation is strong and increases with frequency. As a preliminary study, we compare synthetic wavelet envelopes with the average of FD simulation wavelet envelopes in 50 synthesized random media, which are characterized by the RMS fractional velocity fluctuation = 0.05, correlation scale a  = 5 km and the background wave velocity V 0 = 4 km s –1 . We use the radiation of a 2 Hz Ricker wavelet from a point source. For all the cases of von Kármán order  = 0.1, 0.5 and 1, we find the synthetic wavelet envelopes are a good match to the characteristics of FD simulation wavelet envelopes in a time window starting from the onset through the maximum peak to the time when the amplitude decreases to half the peak amplitude.

Description: We present a thermodynamically based formulation for modelling dynamic rupture processes in the brittle crust using a continuum damage-breakage rheology. The model combines aspects of a continuum viscoelastic damage framework for brittle solids with a continuum breakage mechanics for granular flow within dynamically generated slip zones. The formulation accounts for the density of distributed cracking and other internal flaws in damaged rocks with a scalar damage parameter, and addresses the grain size distribution of a granular phase in the slip zone with a breakage parameter. A dynamic brittle instability is associated with a critical level of damage in the solid, leading to loss of convexity of the solid strain energy, localization and transition to a granular phase associated with lower energy level. The continuum damage-breakage rheology model treats the localization to a slip zone at the onset of dynamic rupture and post-failure recovery process as phase transitions between solid and granular states. The model generates sub- and supershear rupture velocities and pulse-type ruptures seen also in frictional models, and additional important features such as strong dynamic changes of volumetric strain near the rupture front and diversity of nucleation mechanisms. The propagation of rupture front and slip accumulation at a point are correlated with sharp dynamic dilation followed by a gradual decay to a level associated with the final volumetric change associated with the granular phase transition in the slipping zone. The local brittle failure process associated with the solid–granular transition is expected to produce isotropic radiation in addition to the deviatoric terms. The framework significantly extends the ability to model brittle processes in complex geometrical structures and allows analysing the roles of gouge thickness and other parameters on nucleation, rupture and radiation characteristics.

Description: Numerical modeling of seismic waves in heterogeneous porous reservoir rocks is an important tool for the interpretation of seismic surveys in reservoir engineering. We apply globally optimal implicit staggered-grid finite differences (FD) to model 2-D wave propagation in heterogeneous poroelastic media at a low-frequency range (〈10 kHz). We validate the numerical solution by comparing it to an analytical-transient solution obtaining clear seismic wavefields including fast P and slow P and S waves (for a porous media saturated with fluid). The numerical dispersion and stability conditions are derived using von Neumann analysis, showing that over a wide range of porous materials the Courant condition governs the stability and this optimal implicit scheme improves the stability of explicit schemes. High-order explicit FD can be replaced by some lower order optimal implicit FD so computational cost will not be as expensive while maintaining the accuracy. Here, we compute weights for the optimal implicit FD scheme to attain an accuracy of = 10 –8 . The implicit spatial differentiation involves solving tridiagonal linear systems of equations through Thomas’ algorithm.

Description: S -wave receiver functions (SRFs) are used to investigate crustal and upper-mantle structure beneath several ice-covered areas of Antarctica. Moho S-to-P (Sp) arrivals are observed at ~6–8 s in SRF stacks for stations in the Gamburtsev Mountains (GAM) and Vostok Highlands (VHIG), ~5–6 s for stations in the Transantarctic Mountains (TAM) and the Wilkes Basin (WILK), and ~3–4 s for stations in the West Antarctic Rift System (WARS) and the Marie Byrd Land Dome (MBLD). A grid search is used to model the Moho Sp conversion time with Rayleigh wave phase velocities from 18 to 30 s period to estimate crustal thickness and mean crustal shear wave velocity. The Moho depths obtained are between 43 and 58 km for GAM, 36 and 47 km for VHIG, 39 and 46 km for WILK, 39 and 45 km for TAM, 19 and 29 km for WARS and 20 and 35 km for MBLD. SRF stacks for GAM, VHIG, WILK and TAM show little evidence of Sp arrivals coming from upper-mantle depths. SRF stacks for WARS and MBLD show Sp energy arriving from upper-mantle depths but arrival amplitudes do not rise above bootstrapped uncertainty bounds. The age and thickness of the crust is used as a heat flow proxy through comparison with other similar terrains where heat flow has been measured. Crustal structure in GAM, VHIG and WILK is similar to Precambrian terrains in other continents where heat flow ranges from ~41 to 58 mW m –2 , suggesting that heat flow across those areas of East Antarctica is not elevated. For the WARS, we use the Cretaceous Newfoundland–Iberia rifted margins and the Mesozoic-Tertiary North Sea rift as tectonic analogues. The low-to-moderate heat flow reported for the Newfoundland–Iberia margins (40–65 mW m –2 ) and North Sea rift (60–85 mW m –2 ) suggest that heat flow across the WARS also may not be elevated. However, the possibility of high heat flow associated with localized Cenozoic extension or Cenozoic-recent magmatic activity in some parts of the WARS cannot be ruled out.

Description: The 2014 April 1, M w 8.1 Iquique (Chile) earthquake struck in the Northern Chile seismic gap. With a rupture length of less than 200 km, it left unbroken large segments of the former gap. Early studies were able to model the main rupture features but results are ambiguous with respect to the role of aseismic slip and left open questions on the remaining hazard at the Northern Chile gap. A striking observation of the 2014 earthquake has been its extensive preparation phase, with more than 1300 events with magnitude above M L 3, occurring during the 15 months preceding the main shock. Increasing seismicity rates and observed peak magnitudes accompanied the last three weeks before the main shock. Thanks to the large data sets of regional recordings, we assess the precursor activity, compare foreshocks and aftershocks and model rupture preparation and rupture effects. To tackle inversion challenges for moderate events with an asymmetric network geometry, we use full waveforms techniques to locate events, map the seismicity rate and derive source parameters, obtaining moment tensors for more than 300 events (magnitudes M w 4.0–8.1) in the period 2013 January 1–2014 April 30. This unique data set of fore- and aftershocks is investigated to distinguish rupture process models and models of strain and stress rotation during an earthquake. Results indicate that the spatial distributions of foreshocks delineated the shallower part of the rupture areas of the main shock and its largest aftershock, well matching the spatial extension of the aftershocks cloud. Most moment tensors correspond to almost pure double couple thrust mechanisms, consistent with the slab orientation. Whereas no significant differences are observed among thrust mechanisms in different areas, nor among thrust foreshocks and aftershocks, the early aftershock sequence is characterized by the presence of normal fault mechanisms, striking parallel to the trench but dipping westward. These events likely occurred in the shallow wedge structure close to the slab interface and are consequence of the increased extensional stress in this region after the largest events. The overall stress inversion result suggests a minor stress rotation after the main shock, but a significant release of the deviatoric stress. The temporal change in the distribution of focal mechanisms can also be explained in terms of the spatial heterogeneity of the stress field: under such interpretation, the potential of a large megathrust earthquake breaking a larger segment offshore Northern Chile remains high.

Description: We analyse the performance of a higher order accurate staggered viscoelastic time-domain finite-difference method, in which the staggered Adams–Bashforth (ABS) third-order and fourth-order accurate time integrators are used for temporal discretization. ABS is a multistep method that uses previously calculated wavefields to increase the order of accuracy in time. The analysis shows that the numerical dispersion is much lower than that of the widely used second-order leapfrog method. Numerical dissipation is introduced by the ABS method which is significantly smaller for fourth-order than third-order accuracy. In 1-D and 3-D simulation experiments, we verify the convincing improvements of simulation accuracy of the fourth-order ABS method. In a realistic elastic 3-D scenario, the computing time reduces by a factor of approximately 2.4, whereas the memory requirements increase by approximately a factor of 2.2. The ABS method thus provides an alternative strategy to increase the simulation accuracy in time by investing computer memory instead of computing time.

Description: Synthetic dispersion curves are here computed in the frame of an ambient-vibration full-wavefield model, which relies on the description of both ambient-vibration ground displacement and its sources as stochastic fields defined on the Earth's surface, stationary in time and homogeneous in space. In this model, previously developed for computing synthetic Horizontal-to-Vertical Spectral Ratio curves, the power spectral density function and the spatial autocorrelation of the displacement are naturally described as functions of the power spectral density function of the generating forces and of the subsoil properties ( via the relevant Green's function), by also accounting for spatial correlation of these forces. Dispersion curves are computed from the displacement power spectral density function and from the spatial autocorrelation according with the well-known f-k and SPAC techniques, respectively. Two examples illustrate the way this new ambient-vibration model works, showing its possible use in better understanding the role of the surface waves in forming the dispersion curves, as well as its capability to capture some remarkable experimental findings.

Description: The Abruzzi region is located in the Central Apennines Neogene fold-and-thrust belt and has one of the highest seismogenic potential in Italy, with high and diffuse crustal seismicity related to NE–SW oriented extension. In this study, we investigate the detailed spatial variation in shear wave splitting providing high-resolution anisotropic structure beneath the L'Aquila region. To accomplish this, we performed a systematic analysis of crustal anisotropic parameters: fast polarization direction ( ) and delay time ( t ). We benefit from the dense coverage of seismic stations operating in the area and from a catalogue of several accurate earthquake locations of the 2009 L'Aquila seismic sequence, related to the M w 6.1 2009 L'Aquila main shock, to describe in detail the geometry of the anisotropic volume around the active faults that ruptured. The spatial variations both in and t suggest a complex anisotropic structure beneath the region caused by a combination of both structural- and stress-induced mechanisms. The average is NNW–SSE oriented (N141°), showing clear similarity both with the local fault strike and the SH max . In the central part of the study area fast axes are oriented NW–SE, while moving towards the northeastern and northwestern sectors the fast directions clearly diverge from the general trend of NW–SE and rotate accordingly to the local fault strikes. The above-mentioned fault-parallel distribution suggests that the observed anisotropy is mostly controlled by the local fault-related structure. Toward the southeast fast directions become orthogonal both to strike of the local mapped faults and to the SH max . Here, are predominantly oriented NE–SW; we interpret this orientation as due to the presence of a highly fractured and overpressurized rock volume which should be responsible of the 90° flips in and the increase in t . Another possible mechanism for NE–SW orientation of in the southeastern sector could be ascribed to the presence of a buried, deep NE–SW oriented fault system. t , both unnormalized and normalized, does not show any clear evidence of increasing with increasing depth indicating that the anisotropy is confined primarily to the shallower crustal layers (~10 km depth). Interpolating t show that higher values are found at the edges of the main patches of the rupture related to the 2009 main shock, while lower values are limited in the central part of the fault plane, where the coseismic slip was higher. We infer that in the areas surrounding the ruptured region, lateral variations in material properties caused overpressurized fluid conditions, while within the main shock ruptured area, high energy released produced an open crack system such that overpressurization was not possible.

Description: This paper presents a set of deterministic 3-D ground motion simulations for the greater metropolitan area of Medellín in the Aburrá Valley, an earthquake-prone region of the Colombian Andes that exhibits moderate-to-strong topographic irregularities. We created the velocity model of the Aburrá Valley region (version 1) using the geological structures as a basis for determining the shear wave velocity. The irregular surficial topography is considered by means of a fictitious domain strategy. The simulations cover a 50 x 50 x 25 km 3 volume, and four M w = 5 rupture scenarios along a segment of the Romeral fault, a significant source of seismic activity in Colombia. In order to examine the sensitivity of ground motion to the irregular topography and the 3-D effects of the valley, each earthquake scenario was simulated with three different models: (i) realistic 3-D velocity structure plus realistic topography, (ii) realistic 3-D velocity structure without topography, and (iii) homogeneous half-space with realistic topography. Our results show how surface topography affects the ground response. In particular, our findings highlight the importance of the combined interaction between source-effects, source-directivity, focusing, soft-soil conditions, and 3-D topography. We provide quantitative evidence of this interaction and show that topographic amplification factors can be as high as 500 per cent at some locations. In other areas within the valley, the topographic effects result in relative reductions, but these lie in the 0–150 per cent range.

Description: Wave-induced local fluid flow is known as a key mechanism to explain the intrinsic wave dissipation in fluid-saturated rocks. Understanding the relationship between the acoustic properties of rocks and fluid patch distributions is important to interpret the observed seismic wave phenomena. A triple-layer patchy (TLP) model is proposed to describe the P -wave dissipation process in a double-porosity media saturated with two immiscible fluids. The double-porosity rock consists of a solid matrix with unique host porosity and inclusions which contain the second type of pores. Two immiscible fluids are considered in concentric spherical patches, where the inner pocket and the outer sphere are saturated with different fluids. The kinetic and dissipation energy functions of local fluid flow (LFF) in the inner pocket are formulated through oscillations in spherical coordinates. The wave propagation equations of the TLP model are based on Biot's theory and the corresponding Lagrangian equations. The P -wave dispersion and attenuation caused by the Biot friction mechanism and the local fluid flow (related to the pore structure and the fluid distribution) are obtained by a plane-wave analysis from the Christoffel equations. Numerical examples and laboratory measurements indicate that P -wave dispersion and attenuation are significantly influenced by the spatial distributions of both, the solid heterogeneity and the fluid saturation distribution. The TLP model is in reasonably good agreement with White's and Johnson's models. However, differences in phase velocity suggest that the heterogeneities associated with double-porosity and dual-fluid distribution should be taken into account when describing the P -wave dispersion and attenuation in partially saturated rocks.

Description: Seismic ambient noise tomography is applied to central and southern Mozambique, located in the tip of the East African Rift (EAR). The deployment of MOZART seismic network, with a total of 30 broad-band stations continuously recording for 26 months, allowed us to carry out the first tomographic study of the crust under this region, which until now remained largely unexplored at this scale. From cross-correlations extracted from coherent noise we obtained Rayleigh wave group velocity dispersion curves for the period range 5–40 s. These dispersion relations were inverted to produce group velocity maps, and 1-D shear wave velocity profiles at selected points. High group velocities are observed at all periods on the eastern edge of the Kaapvaal and Zimbabwe cratons, in agreement with the findings of previous studies. Further east, a pronounced slow anomaly is observed in central and southern Mozambique, where the rifting between southern Africa and Antarctica created a passive margin in the Mesozoic, and further rifting is currently happening as a result of the southward propagation of the EAR. In this study, we also addressed the question concerning the nature of the crust (continental versus oceanic) in the Mozambique Coastal Plains (MCP), still in debate. Our data do not support previous suggestions that the MCP are floored by oceanic crust since a shallow Moho could not be detected, and we discuss an alternative explanation for its ocean-like magnetic signature. Our velocity maps suggest that the crystalline basement of the Zimbabwe craton may extend further east well into Mozambique underneath the sediment cover, contrary to what is usually assumed, while further south the Kaapval craton passes into slow rifted crust at the Lebombo monocline as expected. The sharp passage from fast crust to slow crust on the northern part of the study area coincides with the seismically active NNE-SSW Urema rift, while further south the Mazenga graben adopts an N-S direction parallel to the eastern limit of the Kaapvaal craton. We conclude that these two extensional structures herald the southward continuation of the EAR, and infer a structural control of the transition between the two types of crust on the ongoing deformation.

Description: Seismic wavefields interact with low-velocity fault damage zones. Waveforms of ballistic fault zone head waves, trapped waves, reflected waves and signatures of trapped noise can provide important information on structural and mechanical fault zone properties. Here we extend the class of observable fault zone waves and reconstruct in-fault reverberations or multiples in a strike-slip faulting environment. Manifestations of the reverberations are significant, consistent wave fronts in the coda of cross-correlation functions that are obtained from scattered earthquake waveforms and seismic noise recorded by a linear fault zone array. The physical reconstruction of Green's functions is evident from the high similarity between the signals obtained from the two different scattered wavefields. Modal partitioning of the reverberation wavefield can be tuned using different data normalization techniques. The results imply that fault zones create their own ambiance, and that the here reconstructed reverberations are a key seismic signature of wear zones. Using synthetic waveform modelling we show that reverberations can be used for the imaging of structural units by estimating the location, extend and magnitude of lateral velocity contrasts. The robust reconstruction of the reverberations from noise records suggests the possibility to resolve the response of the damage zone material to various external and internal loading mechanisms.

Description: The Lower Rhine Embayment in Central Europe hosts a rift system that has very low deformation rates. The faults in this area have slip rates of less than 0.1 mm yr –1 , which does not allow to investigate ongoing tectonic deformation with geodetic techniques, unless they cover very long time spans. Instrumental seismicity does only cover a small fraction of the very long earthquake recurrence intervals of several thousands of years. Palaeoseismological studies are needed to constrain slip rates and the earthquake history of such faults. Destructive earthquakes are rare in the study area, but did occur in historic times. In 1755/1756, a series of strong earthquakes caused significant destruction in the city of Düren (Germany) and the surrounding areas. In this study we document palaeoseismological data from the nearby Rurrand Fault. In contrast to earlier studies on the same fault, we found evidence for a surface rupturing earthquake in the Holocene, and we identified at least one more surface rupturing event. Our study shows that the Rurrand Fault currently accommodates deformation in earthquakes rather than by creeping. The coseismic offsets were determined to be between less than 0.5 m per event. We assign maximum possible magnitudes of M w 5.9–6.8 for the Rurrand Fault and a slip rate of at least 0.02–0.03 mm yr –1 for the last ~130–50 kyr. The surface ruptures did not occur at the main fault trace that has a clear morphological expression due to older tectonic motions, but on a younger fault strand in the hanging wall of the main fault. Terrain analyses based on 1 m resolution airborne LiDAR data have been used to image the subtle morphological expression of this young fault zone. Georadar and electric resistivity tomography were applied to image the fault zone at depth and to test if these shallow geophysical methods can be used to identify and trace the fault zone. Georadar failed to produce reliable results, but geoelectrics were successfully applied and allowed us to retrieve slip rate estimates. Our results indicate that the Düren 1755/1756 earthquakes did not produce surface ruptures at the Rurrand Fault, either because they did not rupture the surface at all, or because they occurred at another, neighbouring fault.

Description: We have measured radial anisotropy in the crust beneath the northern part of North China by jointly inverting Rayleigh and Love wave phase velocities at periods less than 35 s from 14 months of ambient noise data recorded by 222 broad-band seismic stations. We also estimate the azimuthal anisotropy of phase velocity from Rayleigh wave data. The fast direction of azimuthal anisotropy varies with periods, NE–SW orientation at short and intermediate periods (10–16 s) and NW–SE orientation at periods larger than 20 s. The NE–SW oriented fast direction of azimuthal anisotropy may be related to the fossilized structural fabrics due to the compression during the Indosinian orogeny from late Palaeozoic to middle Mesozoic. The NW–SE trend of anisotropic fabric in the lower crust and uppermost mantle is probably associated with the later lithospheric extension. The observed radial anisotropy also shows a two-layer feature, negative radial anisotropy (Vsh 〈 Vsv) in the upper crust and positive (Vsh 〉 Vsv) in the middle-lower crust. The compressional tectonics from late Palaeozoic to middle Mesozoic may cause crustal materials align vertically throughout the crust. This vertical fabric could make Vsh slower than Vsv. However, the lithospheric extension in the late Mesozoic to Cenozoic time could overprint the older fabric in the middle and lower crust by magma intrusion and underplating. Horizontal alignment of the material or intruded melt sills due to the extension probably produce the observed strong positive radial anisotropy in the middle and lower crust.

Description: Using active and passive seismology data we derive a shear ( S ) wave velocity model and a Poisson's ratio () model across the Chilean convergent margin along a profile at 38°15'S, where the M w 9.5 Valdivia earthquake occurred in 1960. The derived S -wave velocity model was constructed using three independently obtained velocity models that were merged together. In the upper part of the profile (0–2 km depth), controlled source data from explosions were used to obtain an S -wave traveltime tomogram. For the middle part (2–20 km depth), data from a temporary seismology array were used to carry out a dispersion analysis. The resulting dispersion curves were used to obtain a 3-D S -wave velocity model. In the lower part (20–75 km depth, depending on the longitude), an already existent local earthquake tomographic image was merged with the other two sections. This final S -wave velocity model and already existent compressional ( P ) wave velocity models along the same transect allowed us to obtain a Poisson's ratio model. The results of this study show that the velocities and Poisson's ratios in the continental crust of this part of the Chilean convergent margin are in agreement with geological features inferred from other studies and can be explained in terms of normal rock types. There is no requirement to call on the existence of measurable amounts of present-day fluids, in terms of seismic velocities, above the plate interface in the continental crust of the Coastal Cordillera and the Central Valley in this part of the Chilean convergent margin. This is in agreement with a recent model of water being transported down and released from the subduction zone.

Description: In the estimation of spatiotemporal slip models, kinematic source inversions using Akaike's Bayesian Information Criterion (ABIC) and the multiple-time-window method have often been used. However, there are cases in which conventional ABIC-based source inversions do not work well in the determination of hyperparameters when a non-negative slip constraint is used. In order to overcome this problem, a new source inversion method was developed in this study. The new method introduces a fully Bayesian method into the kinematic multiple-time-window source inversion. The multiple-time-window method is one common way of parametrizing a source time function and is highly flexible in terms of the shape of the source time function. The probability distributions of model parameters and hyperparameters can be directly obtained by using the Markov chain Monte Carlo method. These probability distributions are useful for simply evaluating the uniqueness and reliability of the derived model, which is another advantage of a fully Bayesian method. This newly developed source inversion method was applied to the 2011 Ibaraki-oki, Japan, earthquake ( M w 7.9) to demonstrate its usefulness. It was demonstrated that the problem with using the conventional ABIC-based source inversion method for hyperparameter determination appeared in the spatiotemporal source inversion of this event and that the newly developed source inversion could overcome this problem.

Description: Disease-gene identification is a challenging process that has multiple applications within functional genomics and personalized medicine. Typically, this process involves both finding genes known to be associated with the disease (through literature search) and carrying out preliminary experiments or screens (e.g. linkage or association studies, copy number analyses, expression profiling) to determine a set of promising candidates for experimental validation. This requires extensive time and monetary resources. We describe Beegle , an online search and discovery engine that attempts to simplify this process by automating the typical approaches. It starts by mining the literature to quickly extract a set of genes known to be linked with a given query, then it integrates the learning methodology of Endeavour (a gene prioritization tool) to train a genomic model and rank a set of candidate genes to generate novel hypotheses. In a realistic evaluation setup, Beegle has an average recall of 84% in the top 100 returned genes as a search engine, which improves the discovery engine by 12.6% in the top 5% prioritized genes. Beegle is publicly available at http://beegle.esat.kuleuven.be/ .

Description: Alternative splicing is an important mechanism in eukaryotes that expands the transcriptome and proteome significantly. It plays an important role in a number of biological processes. Understanding its regulation is hence an important challenge. Recently, increasing evidence has been collected that supports an involvement of intragenic DNA methylation in the regulation of alternative splicing. The exact mechanisms of regulation, however, are largely unknown, and speculated to be complex: different methylation profiles might exist, each of which could be associated with a different regulation mechanism. We present a computational technique that is able to determine such stable methylation patterns and allows to correlate these patterns with inclusion propensity of exons. Pattern detection is based on dynamic time warping (DTW) of methylation profiles, a sophisticated similarity measure for signals that can be non-trivially transformed. We design a flexible self-organizing map approach to pattern grouping. Exemplary application on available data sets indicates that stable patterns which correlate non-trivially with exon inclusion do indeed exist. To improve the reliability of these predictions, further studies on larger data sets will be required. We have thus taken great care that our software runs efficiently on modern hardware, so that it can support future studies on large-scale data sets.

Description: In a set of problems ranging from 4-D seismic to salt boundary estimation, updates to the velocity model often have a highly localized nature. Numerical techniques for these applications such as full-waveform inversion (FWI) require an estimate of the wavefield to compute the model updates. When dealing with localized problems, it is wasteful to compute these updates in the global domain, when we only need them in our region of interest. This paper introduces a local solver that generates forward and adjoint wavefields which are, to machine precision, identical to those generated by a full-domain solver evaluated within the region of interest. This means that the local solver computes all interactions between model updates within the region of interest and the inhomogeneities in the background model outside. Because no approximations are made in the calculation of the forward and adjoint wavefields, the local solver can compute the identical gradient in the region of interest as would be computed by the more expensive full-domain solver. In this paper, the local solver is used to efficiently generate the FWI gradient at the boundary of a salt body. This gradient is then used in a level set method to automatically update the salt boundary.

Description: Full waveform inversion using the conventional L 2 distance to measure the misfit between seismograms is known to suffer from cycle skipping. An alternative strategy is proposed in this study, based on a measure of the misfit computed with an optimal transport distance. This measure allows to account for the lateral coherency of events within the seismograms, instead of considering each seismic trace independently, as is done generally in full waveform inversion. The computation of this optimal transport distance relies on a particular mathematical formulation allowing for the non-conservation of the total energy between seismograms. The numerical solution of the optimal transport problem is performed using proximal splitting techniques. Three synthetic case studies are investigated using this strategy: the Marmousi 2 model, the BP 2004 salt model, and the Chevron 2014 benchmark data. The results emphasize interesting properties of the optimal transport distance. The associated misfit function is less prone to cycle skipping. A workflow is designed to reconstruct accurately the salt structures in the BP 2004 model, starting from an initial model containing no information about these structures. A high-resolution P -wave velocity estimation is built from the Chevron 2014 benchmark data, following a frequency continuation strategy. This estimation explains accurately the data. Using the same workflow, full waveform inversion based on the L 2 distance converges towards a local minimum. These results yield encouraging perspectives regarding the use of the optimal transport distance for full waveform inversion: the sensitivity to the accuracy of the initial model is reduced, the reconstruction of complex salt structure is made possible, the method is robust to noise, and the interpretation of seismic data dominated by reflections is enhanced.

Description: This paper concerns the numerical modelling of time-domain mechanical waves in viscoelastic media based on a generalized Zener model. To do so, classically in the literature relaxation mechanisms are introduced, resulting in a set of the so-called memory variables and thus in large computational arrays that need to be stored. A challenge is thus to accurately mimic a given attenuation law using a minimal set of relaxation mechanisms. For this purpose, we replace the classical linear approach of Emmerich & Korn with a nonlinear optimization approach with constraints of positivity. We show that this technique is more accurate than the linear approach. Moreover, it ensures that physically meaningful relaxation times that always honour the constraint of decay of total energy with time are obtained. As a result, these relaxation times can always be used in a stable way in a modelling algorithm, even in the case of very strong attenuation for which the classical linear approach may provide some negative and thus unusable coefficients.

Description: Fractures are a source of extra compliance in the rock mass. The mechanical properties of the rock matrix as well as the propagation of seismic waves inside the rock medium are dependent on the magnitude of roughness and offset between the imperfect fracture interfaces. Fracture compliance can estimate the degree of contact between fracture faces, type of fluid filling the fracture and the fracture roughness. To characterize these fracture properties, compliance ratio, known by the ratio of normal-to-shear compliance, can be a potential tool in the subsurface studies to improve the well layout design. The focus of this study is to illustrate how the compliance ratio of a rough fracture, with or without the offset between the fracture faces, can diverge from the compliance ratio of a fracture with smooth interface. Quasi-static and dynamic methods are two common ways to calculate the compliance. The former calculates the compliance by measuring the change in the displacement with the applied stress, while the latter estimates the compliance through monitoring the changes in propagation of seismic waves. To compare the compliance ratios of fractures with imperfect and smooth interfaces in an infinite medium, a numerical finite-element model is built in commercial finite-element software. The imperfect interface of the fracture is modeled with saw-tooth-like structures where they can be partially or fully in contact. The defined saw-tooth-like structures of contact asperities impose an in-plane asymmetry in the shear direction. This asymmetry causes two different values for the compliance in shear direction, known as the soft and stiff shear compliance. Our numerical simulations suggest the increase in the degree of contact between the fracture faces increases the compliance ratio in the stiff direction more than the soft direction. The compliance ratio of the fracture with the imperfect interface is larger than the compliance ratio of the smooth fracture. We suggest that the interlocking and riding up effects at the fracture interface may explain our findings in this study.

Description: Ambient noise correlation method has been widely used to extract surface waves and tomography. The extraction of body waves has been very limited, but recent reports have suggested promises for deep incident waves. Here we report our first observations of triplicated PKP phases (important phases for studying the Earth's core) and confirm observations of other body-wave core phases from noise correlations. We use dense seismic arrays in South America and China Regional Seismic Networks at distances from 145° to the antipode. We can clearly observe different PKP branches (df, bc and ab) in stacks of the station–station correlations. Both ambient noise and earthquake coda contribute to PKP phases. However, the contributions vary with frequency and with body-wave phases. At shorter periods (5–20 s), three branches of PKP (df, bc and ab) can be extracted from ambient noise and the ab phase from earthquake coda. At longer periods (15–50 s), earthquake coda are effective in generating the df branch, but not the ab branch. The generation of the PKIKP phase (df branch) from earthquake coda does not depend on earthquake focal mechanisms or focal depths. However, earthquakes far from the stations contribute more than events closer by. The best coda window is around 10 000–40 000 s and the best magnitude threshold is M w greater than 6.8 or 6.9. The observation of triplicated PKP branches from noise correlations provides a new type of data for studying the Earth's deep interior, in particularly the inner core anisotropy, which overcomes some of the limitations of traditional earthquake-based studies (such as limited source distributions and source location errors).

Description: The project Seismic Hazard Harmonization in Europe (SHARE), completed in 2013, presents significant improvements over previous regional seismic hazard modeling efforts. The Global Strain Rate Map v2.1, sponsored by the Global Earthquake Model Foundation and built on a large set of self-consistent geodetic GPS velocities, was released in 2014. To check the SHARE seismic source models that were based mainly on historical earthquakes and active fault data, we first evaluate the SHARE historical earthquake catalogues and demonstrate that the earthquake magnitudes are acceptable. Then, we construct an earthquake potential model using the Global Strain Rate Map data. SHARE models provided parameters from which magnitude–frequency distributions can be specified for each of 437 seismic source zones covering most of Europe. Because we are interested in proposed magnitude limits, and the original zones had insufficient data for accurate estimates, we combine zones into five groups according to SHARE's estimates of maximum magnitude. Using the strain rates, we calculate tectonic moment rates for each group. Next, we infer seismicity rates from the tectonic moment rates and compare them with historical and SHARE seismicity rates. For two of the groups, the tectonic moment rates are higher than the seismic moment rates of the SHARE models. Consequently, the rates of large earthquakes forecast by the SHARE models are lower than those inferred from tectonic moment rate. In fact, the SHARE models forecast higher seismicity rates than the historical rates, which indicate that the authors of SHARE were aware of the potentially higher seismic activities in the zones. For one group, the tectonic moment rate is lower than the seismic moment rates forecast by the SHARE models. As a result, the rates of large earthquakes in that group forecast by the SHARE model are higher than those inferred from tectonic moment rate, but lower than what the historical data show. For the other two groups, the seismicity rates from tectonic moment rate, historical data and SHARE models are consistent. For four groups, the maximum magnitudes used by SHARE are fairly consistent with the probable maximum magnitudes inferred from tectonic strain rates. This study demonstrates that: (1) the strain-rate data are useful for constraining seismicity rates and magnitude limits; and (2) SHARE seismic source models and historical earthquake catalogues are satisfactory.

Description: Tumors are characterized by properties of genetic instability, heterogeneity, and significant oligoclonality. Elucidating this intratumoral heterogeneity is challenging but important. In this study, we propose a framework, BubbleTree, to characterize the tumor clonality using next generation sequencing (NGS) data. BubbleTree simultaneously elucidates the complexity of a tumor biopsy, estimating cancerous cell purity, tumor ploidy, allele-specific copy number, and clonality and represents this in an intuitive graph. We further developed a three-step heuristic method to automate the interpretation of the BubbleTree graph, using a divide-and-conquer strategy. In this study, we demonstrated the performance of BubbleTree with comparisons to similar commonly used tools such as THetA2, ABSOLUTE, AbsCN-seq and ASCAT, using both simulated and patient-derived data. BubbleTree outperformed these tools, particularly in identifying tumor subclonal populations and polyploidy. We further demonstrated BubbleTree's utility in tracking clonality changes from patients’ primary to metastatic tumor and dating somatic single nucleotide and copy number variants along the tumor clonal evolution. Overall, the BubbleTree graph and corresponding model is a powerful approach to provide a comprehensive spectrum of the heterogeneous tumor karyotype in human tumors. BubbleTree is R-based and freely available to the research community ( https://www.bioconductor.org/packages/release/bioc/html/BubbleTree.html ).

Description: We present a method for estimating ratios of P and S waves corner frequencies ( R cf ) and earthquake strain drops by joint analysis of P and S source spectra of neighbouring groups of events. The method is applied systematically to data generated by ~9000 earthquakes around the Karadere segment of the North Anatolian Fault Zone. The results indicate several regions that produce consistently R cf values higher (e.g. 〉2) than expected from classical earthquake source models. These are associated generally with fault sections having strong geometrical heterogeneities, shallow depth sections and/or locations without large pre-existing surface trace. Earthquake ruptures in such regions are likely to generate significant rock damage and tensile components of faulting. To assess whether the observed high R cf values are produced by enriched high frequency P waves, reduced high frequency S waves or both, we compare the associated P and S spectra with mean/median results. The analysis suggests that the high R cf values of shallow events (depth 〈4 km) are generated primarily by reduced high frequency S radiation, and that the contribution from elevated high frequency P radiation increases with depth and proximity to geometrical complexities. The results highlight the importance of considering carefully the existence of some volumetric source components in earthquake rupture processes.

Description: Fractures are a source of extra compliance in the rock mass. The mechanical properties of the rock matrix as well as the propagation of seismic waves inside the rock medium are dependent on the magnitude of roughness and offset between the imperfect fracture interfaces. Fracture compliance can estimate the degree of contact between fracture faces, type of fluid filling the fracture and the fracture roughness. To characterize these fracture properties, compliance ratio, known by the ratio of normal-to-shear compliance, can be a potential tool in the subsurface studies to improve the well layout design. The focus of this study is to illustrate how the compliance ratio of a rough fracture, with or without the offset between the fracture faces, can diverge from the compliance ratio of a fracture with smooth interface. Quasi-static and dynamic methods are two common ways to calculate the compliance. The former calculates the compliance by measuring the change in the displacement with the applied stress, while the latter estimates the compliance through monitoring the changes in propagation of seismic waves. To compare the compliance ratios of fractures with imperfect and smooth interfaces in an infinite medium, a numerical finite-element model is built in commercial finite-element software. The imperfect interface of the fracture is modeled with saw-tooth-like structures where they can be partially or fully in contact. The defined saw-tooth-like structures of contact asperities impose an in-plane asymmetry in the shear direction. This asymmetry causes two different values for the compliance in shear direction, known as the soft and stiff shear compliance. Our numerical simulations suggest the increase in the degree of contact between the fracture faces increases the compliance ratio in the stiff direction more than the soft direction. The compliance ratio of the fracture with the imperfect interface is larger than the compliance ratio of the smooth fracture. We suggest that the interlocking and riding up effects at the fracture interface may explain our findings in this study.

Description: Full waveform inversion using the conventional L 2 distance to measure the misfit between seismograms is known to suffer from cycle skipping. An alternative strategy is proposed in this study, based on a measure of the misfit computed with an optimal transport distance. This measure allows to account for the lateral coherency of events within the seismograms, instead of considering each seismic trace independently, as is done generally in full waveform inversion. The computation of this optimal transport distance relies on a particular mathematical formulation allowing for the non-conservation of the total energy between seismograms. The numerical solution of the optimal transport problem is performed using proximal splitting techniques. Three synthetic case studies are investigated using this strategy: the Marmousi 2 model, the BP 2004 salt model, and the Chevron 2014 benchmark data. The results emphasize interesting properties of the optimal transport distance. The associated misfit function is less prone to cycle skipping. A workflow is designed to reconstruct accurately the salt structures in the BP 2004 model, starting from an initial model containing no information about these structures. A high-resolution P -wave velocity estimation is built from the Chevron 2014 benchmark data, following a frequency continuation strategy. This estimation explains accurately the data. Using the same workflow, full waveform inversion based on the L 2 distance converges towards a local minimum. These results yield encouraging perspectives regarding the use of the optimal transport distance for full waveform inversion: the sensitivity to the accuracy of the initial model is reduced, the reconstruction of complex salt structure is made possible, the method is robust to noise, and the interpretation of seismic data dominated by reflections is enhanced.

Description: We present a novel earthquake location method using acoustic wave-equation-based traveltime inversion. The linear relationship between the location perturbation ( t 0 , x s ) and the resulting traveltime residual t of a particular seismic phase, represented by the traveltime sensitivity kernel K ( t 0 , x s ) with respect to the earthquake location ( t 0 , x s ), is theoretically derived based on the adjoint method. Traveltime sensitivity kernel K ( t 0 , x s ) is formulated as a convolution between the forward and adjoint wavefields, which are calculated by numerically solving two acoustic wave equations. The advantage of this newly derived traveltime kernel is that it not only takes into account the earthquake–receiver geometry but also accurately honours the complexity of the velocity model. The earthquake location is obtained by solving a regularized least-squares problem. In 3-D realistic applications, it is computationally expensive to conduct full wave simulations. Therefore, we propose a 2.5-D approach which assumes the forward and adjoint wave simulations within a 2-D vertical plane passing through the earthquake and receiver. Various synthetic examples show the accuracy of this acoustic wave-equation-based earthquake location method. The accuracy and efficiency of the 2.5-D approach for 3-D earthquake location are further verified by its application to the 2004 Big Bear earthquake in Southern California.

Description: Seismic interferometry applied to 120 hr of railroad traffic recorded by an array of vertical component seismographs along a railway within the Rio Grande rift has recovered surface and body waves characteristic of the geology beneath the railway. Linear and hyperbolic arrivals are retrieved that agree with surface (Rayleigh), direct and reflected P waves observed by nearby conventional seismic surveys. Train-generated Rayleigh waves span a range of frequencies significantly higher than those recovered from typical ambient noise interferometry studies. Direct P -wave arrivals have apparent velocities appropriate for the shallow geology of the survey area. Significant reflected P -wave energy is also present at relatively large offsets. A common midpoint stack produces a reflection image consistent with nearby conventional reflection data. We suggest that for sources at the free surface (e.g. trains) increasing the aperture of the array to record wide angle reflections, in addition to longer recording intervals, might allow the recovery of deeper geological structure from railroad traffic. Frequency–wavenumber analyses of these recordings indicate that the train source is symmetrical (i.e. approaching and receding) and that deeper refracted energy is present although not evident in the time-offset domain. These results confirm that train-generated vibrations represent a practical source of high-resolution subsurface information, with particular relevance to geotechnical and environmental applications.

Description: An approach to depth migration, based on an integral representation of seismic data, that is, wavefields recorded on the boundary, is presented in terms of Poincaré wavelets. Each wavelet is taken as a boundary datum for a high-frequency asymptotic solution of the wave equation. This solution, which we call the quasiphoton or the Gaussian wave packet, decreases in a Gaussian manner away from a point running along a ray that is launched from the surface. The deformation of the propagating packet is taken into account in the migration algorithm. A numerical example of zero-offset migration with synthetic seismograms calculated for the 2-D SEG/EAGE salt model is presented. The result, which uses only 3.9 per cent of the total number of coefficients, is a satisfactory image, with a threshold of 0.75 per cent.

Description: The motion of a self-gravitating hyperelastic body is described through a time-dependent mapping from a reference body into physical space, and its material properties are determined by a referential density and strain-energy function defined relative to the reference body. Points within the reference body do not have a direct physical meaning, but instead act as particle labels that could be assigned in different ways. We use Hamilton's principle to determine how the referential density and strain-energy functions transform when the particle labels are changed, and describe an associated ‘particle relabelling symmetry’. We apply these results to linearized elastic wave propagation and discuss their implications for seismological inverse problems. In particular, we show that the effects of boundary topography on elastic wave propagation can be mapped exactly into volumetric heterogeneity while preserving the form of the equations of motion. Several numerical calculations are presented to illustrate our results.

Description: Scattering due to small-scale heterogeneities in the rock formation surrounding a wellbore can significantly change the acoustic waveform from a logging measurement which in turn can be used to characterize the formation heterogeneities. This study simulates the elastic heterogeneity scattering in monopole and dipole acoustic logging and analyse the resulting effects on the waveforms. The results show that significant coda waves are generated in both monopole and dipole waveforms and the dipole coda is dominated by S -to- S scattering, which can be effectively utilized to diagnose the heterogeneity in the rock formation. The coda wave modelling and analysis were used to characterize dipole acoustic data logged before and after fracturing a reservoir interval, with significant coda wave in the after-fracturing data indicating fracturing-induced heterogeneous property change in the rock volume surrounding the borehole.

Description: We present a novel earthquake location method using acoustic wave-equation-based traveltime inversion. The linear relationship between the location perturbation ( t 0 , x s ) and the resulting traveltime residual t of a particular seismic phase, represented by the traveltime sensitivity kernel K ( t 0 , x s ) with respect to the earthquake location ( t 0 , x s ), is theoretically derived based on the adjoint method. Traveltime sensitivity kernel K ( t 0 , x s ) is formulated as a convolution between the forward and adjoint wavefields, which are calculated by numerically solving two acoustic wave equations. The advantage of this newly derived traveltime kernel is that it not only takes into account the earthquake–receiver geometry but also accurately honours the complexity of the velocity model. The earthquake location is obtained by solving a regularized least-squares problem. In 3-D realistic applications, it is computationally expensive to conduct full wave simulations. Therefore, we propose a 2.5-D approach which assumes the forward and adjoint wave simulations within a 2-D vertical plane passing through the earthquake and receiver. Various synthetic examples show the accuracy of this acoustic wave-equation-based earthquake location method. The accuracy and efficiency of the 2.5-D approach for 3-D earthquake location are further verified by its application to the 2004 Big Bear earthquake in Southern California.

Description: Observations of shear wave splitting in SKS seismic phase play a key role in the current efforts to understand kinematics and dynamics of mantle flow, but azimuthal anisotropy as a depth-localized phenomenon still is poorly known. Here we analyse stratification of seismic azimuthal anisotropy beneath central and northern Anatolia (a microplate within the Alpine belt) by inverting P -wave receiver functions jointly with shear wave splitting in SKS seismic phase. The analysis is based on recordings of stations of the North Anatolian Fault (NAF) passive seismic experiment. In the resulting model in a depth interval from 120 to 200 km fast direction of anisotropy is nearly parallel to the plate motion direction (~E–W), whilst a normal direction (close to S–N) is found in the low velocity zone (LVZ) between 60 and 90 km. Our preferred interpretation of these data suggests that the flow in upper mantle is nearly parallel to the Anatolian plate motion direction in the depth range from the LAB to 200 km, but in part of the LVZ fast direction of anisotropy is normal to the direction of shear in the mantle. This relation between anisotropy and shear is known from laboratory experiments with peridotite-type rock containing melt. A similar relation between anisotropy and flow in the LVZ is found in Fennoscandia. These findings may have far-reaching implications for interpreting mantle anisotropy elsewhere.

Description: On 2013 April 13, an inland earthquake of M w 5.8 occurred in Awaji Island, which forms the western boundary of the Osaka sedimentary basin in western Japan. The strong ground motion data were collected from more than 100 stations within the basin and it was found that in the Osaka Plain, the pseudo velocity response spectra at a period of around 6.5 s were significantly larger than at other stations of similar epicentral distance outside the basin. The ground motion lasted longer than 3 min in the Osaka Plain where its bedrock depth spatially varies from approximately 1 to 2 km. We modelled long-period (higher than 2 s) ground motions excited by this earthquake, using the finite difference method assuming a point source, to validate the present velocity structure model and to obtain better constraint of the attenuation factor of the sedimentary part of the basin. The effect of attenuation in the simulation was included in the form of Q ( f ) =  Q 0 ( f / f 0 ), where Q 0 at a reference frequency f 0 was given by a function of the S -wave velocity, Q 0 = α V S . We searched for appropriate Q 0 values by changing α for a fixed value of f 0 = 0.2 Hz. It was found that values of α from 0.2 to 0.5 fitted the observations reasonably well, but that the value of α  = 0.3 performed best. Good agreement between the observed and simulated velocity waveforms was obtained for most stations within the Osaka Basin in terms of both amplitude and ground motion duration. However, underestimation of the pseudo velocity response spectra in the period range of 5–7 s was recognized in the central part of the Osaka Plain, which was caused by the inadequate modelling of later phases or wave packets in this period range observed approximately 2 min after the direct S -wave arrival. We analysed this observed later phase and concluded that it was a Love wave originating from the direction of the east coast of Awaji Island.

Description: Geodetic slip inversions for three major ( M w 〉 7) strike-slip earthquakes (1992 Landers, 1999 Hector Mine and 2010 El Mayor–Cucapah) show a 15–60 per cent reduction in slip near the surface (depth 〈 2 km) relative to the slip at deeper depths (4–6 km). This significant difference between surface coseismic slip and slip at depth has been termed the shallow slip deficit (SSD). The large magnitude of this deficit has been an enigma since it cannot be explained by shallow creep during the interseismic period or by triggered slip from nearby earthquakes. One potential explanation for the SSD is that the previous geodetic inversions lack data coverage close to surface rupture such that the shallow portions of the slip models are poorly resolved and generally underestimated. In this study, we improve the static coseismic slip inversion for these three earthquakes, especially at shallow depths, by: (1) including data capturing the near-fault deformation from optical imagery and SAR azimuth offsets; (2) refining the interferometric synthetic aperture radar processing with non-boxcar phase filtering, model-dependent range corrections, more complete phase unwrapping by SNAPHU (Statistical Non-linear Approach for Phase Unwrapping) assuming a maximum discontinuity and an on-fault correlation mask; (3) using more detailed, geologically constrained fault geometries and (4) incorporating additional campaign global positioning system (GPS) data. The refined slip models result in much smaller SSDs of 3–19 per cent. We suspect that the remaining minor SSD for these earthquakes likely reflects a combination of our elastic model's inability to fully account for near-surface deformation, which will render our estimates of shallow slip minima, and potentially small amounts of interseismic fault creep or triggered slip, which could ‘make up’ a small percentages of the coseismic SSD during the interseismic period. Our results indicate that it is imperative that slip inversions include accurate measurements of near-fault surface deformation to reliably constrain spatial patterns of slip during major strike-slip earthquakes.

Description: We present a mathematically consistent numerical algorithm for the simulation of earthquake rupture with rate-and-state friction. Its main features are adaptive time stepping, a novel algebraic solution algorithm involving nonlinear multigrid and a fixed point iteration for the rate-and-state decoupling. The algorithm is applied to a laboratory scale subduction zone which allows us to compare our simulations with experimental results. Using physical parameters from the experiment, we find a good fit of recurrence time of slip events as well as their rupture width and peak slip. Computations in 3-D confirm efficiency and robustness of our algorithm.

Description: Seismic moment tensors can provide information on the size and orientation of fractures producing acoustic emissions (AEs) and on the stress conditions in the sample. The moment tensor inversion of AEs is, however, a demanding procedure requiring carefully calibrated sensors and accurate knowledge of the velocity model. In field observations, the velocity model is usually isotropic and time independent. In laboratory experiments, the velocity is often anisotropic and time dependent and attenuation might be significant due to opening or closure of microcracks in the sample during loading. In this paper, we study the sensitivity of the moment tensor inversion to anisotropy of P -wave velocities and attenuation. We show that retrieved moment tensors critically depend on anisotropy and attenuation and their neglect can lead to misinterpretations of the source mechanisms. The accuracy of the inversion also depends on the fracturing mode of AEs: tensile events are more sensitive to P -wave anisotropy and attenuation than shear events. We show that geometry of faulting in anisotropic rocks should be studied using the source tensors, since the P - and T -axes of the moment tensors are affected by velocity anisotropy and deviate from the true orientation of faulting. The stronger the anisotropy is, the larger the deviations are. Finally, we prove that the moment tensor inversion applied to a large dataset of AEs can be utilized to provide information on the attenuation parameters of the rock sample. The method is capable of measuring anisotropic attenuation in the sample and allows for detection of dilatant cracking according to the stress regime.

Description: We investigate the lithospheric structure of Cameroon inverting Rayleigh waves obtained from the cross-correlation of ambient seismic noise. We correlate seismic records between 32 broad-band stations and we obtain good quality Rayleigh waves for 310 interstation paths. We measure group velocity dispersion curves from the reconstructed Rayleigh waves in the period range 10–35 s and we invert the group velocities for tomographic images. After the tomography the group velocities are then inverted, together with longer period group velocity measurements from existing literature, to compute a 3-D S -wave velocity model of the Cameroon lithosphere down to 100 km depth. Our results provide an unprecedented mapping of the physical properties of the different crustal units and their correlations with surface geology, as well as with mantle lithospheric variations. The Cameroon Volcanic Line (CVL) appears as a segmented feature exhibiting different physical properties along strike. The active Mt Cameroon volcano is underlain by very low velocities, unlike the other segments of the CVL. The along-strike variations in crustal structure suggest that lateral heterogeneities in lithospheric thickness and physical properties have influenced the location and distribution of magmatism. The crust beneath the Central African Shear Zone exhibits a sizeable low velocity anomaly. The lithosphere beneath Cameroon is characterised by a heterogeneous crust with a relatively constant thickness and a low velocity uppermost mantle at the edge of the Congo Craton. Our results favour processes combining small-scale upwelling at the edge of a thick lithosphere and reactivation of Precambrian basement structures to explain the distribution of Holocene-Recent magmatism and plateau uplift. Our results also indicate that Mt Cameroon and surroundings areas are the most at risk zones for magmatic activity during this stage of CVL development.

Description: Although optimal, computing the moment tensor solution is not always a viable option for the calculation of the size of an earthquake, especially for small events (say, below M w 2.0). Here we show an alternative approach to the calculation of the moment-rate spectra of small earthquakes, and thus of their scalar moments, that uses a network-based calibration of crustal wave propagation. The method works best when applied to a relatively small crustal volume containing both the seismic sources and the recording sites. In this study we present the calibration of the crustal volume monitored by the High-Resolution Seismic Network (HRSN), along the San Andreas Fault (SAF) at Parkfield. After the quantification of the attenuation parameters within the crustal volume under investigation, we proceed to the spectral correction of the observed Fourier amplitude spectra for the 100 largest events in our data set. Multiple estimates of seismic moment for the all events (1811 events total) are obtained by calculating the ratio of rms-averaged spectral quantities based on the peak values of the ground velocity in the time domain, as they are observed in narrowband-filtered time-series. The mathematical operations allowing the described spectral ratios are obtained from Random Vibration Theory (RVT). Due to the optimal conditions of the HRSN, in terms of signal-to-noise ratios, our network-based calibration allows the accurate calculation of seismic moments down to M w 〈 0. However, because the HRSN is equipped only with borehole instruments, we define a frequency-dependent Generalized Free-Surface Effect (GFSE), to be used instead of the usual free-surface constant F = 2. Our spectral corrections at Parkfield need a different GFSE for each side of the SAF, which can be quantified by means of the analysis of synthetic seismograms. The importance of the GFSE of borehole instruments increases for decreasing earthquake's size because for smaller earthquakes the bandwidth available for our calculations is consistently shifted towards higher frequencies.

Description: As a result of the accumulation of high-resolution observation data, 3-D high-fidelity crustal structure data for large domains are becoming available. However, it has been difficult to use such data to perform elastic/viscoelastic crustal deformation analyses in large domains with quality assurance of the numerical simulation that guarantees convergence of the numerical solution with respect to the discretization size because the costs of analysis are significantly high. This paper proposes a method of constructing a high-fidelity crustal structure finite element (FE) model using high-fidelity crustal structure data and fast FE analysis to reduce the costs of analysis (based on automatic FE model generation for parallel computation, OpenMP/MPI hybrid parallel computation on distributed memory computers, a geometric multigrid, variable preconditioning and multiple precision arithmetic). Using the proposed methods, we construct 10 billion degree-of-freedom high-fidelity crustal structure FE models for the entire Japan, and conduct elastic/viscoelastic crustal deformation analysis using this model with enough high accuracy of the numerical simulation.

Description: The EarthScope USArray provides an opportunity to obtain detailed images of the continental upper mantle at an unprecedented scale. The majority of mantle models derived from USArray data to date contain spatial variations in seismic-wave speed; however, in many cases these data sets do not by themselves allow a non-unique interpretation. Joint interpretation of seismic attenuation and velocity models can improve upon the interpretations based only on velocity and provide important constraints on the temperature, composition, melt content, and volatile content of the mantle. The surface wave amplitudes that constrain upper-mantle attenuation are sensitive to factors in addition to attenuation, including the earthquake source excitation, focusing and defocusing by elastic structure, and local site amplification. Because of the difficulty of isolating attenuation from these other factors, little is known about the attenuation structure of the North American upper mantle. In this study, Rayleigh wave traveltime and amplitude in the period range 25–100 s are measured using an interstation cross-correlation technique, which takes advantage of waveform similarity at nearby stations. Several estimates of Rayleigh wave attenuation and site amplification are generated at each period, using different approaches to separate the effects of attenuation and local site amplification on amplitude. It is assumed that focusing and defocusing effects can be described by the Laplacian of the traveltime field. All approaches identify the same large-scale patterns in attenuation, including areas where the attenuation values are likely contaminated by unmodelled focusing and defocusing effects. Regionally averaged attenuation maps are constructed after removal of the contaminated attenuation values, and the variations in intrinsic shear attenuation that are suggested by these Rayleigh wave attenuation maps are explored.

Description: Surface wave magnitude M s for a compilation of 72 nuclear tests detonated in hard rock media for which yields and burial depths have been reported in the literature is shown to scale with yield W as a + b x log[ W ], where a = 2.50 ± 0.08 and b = 0.80 ± 0.05. While the exponent b is consistent with an M s scaling model for fully coupled, normal containment-depth explosions, the intercept a is offset 0.45 magnitude units lower than the model. The cause of offset is important to understand in terms of the explosion source. Hard rock explosions conducted in extensional and compressional stress regimes show similar offsets, an indication that the tectonic setting in which an explosion occurs plays no role causing the offset. The scaling model accounts for the effects of source medium material properties on the generation of 20-s period Rayleigh wave amplitudes. Aided by thorough characterizations of the explosion and tectonic release sources, an extensive analysis of the 1963 October 26 Shoal nuclear test detonated in granite 27 miles southeast of Fallon NV shows that the offset is consistent with the predictions of a material damage source model related to non-linear stress wave interactions with the free surface. This source emits Rayleigh waves with polarity opposite to waves emitted by the explosion. The Shoal results were extended to analyse surface waves from the 1962 February 15 Hardhat nuclear test, the 1988 September 14 Soviet Joint Verification Experiment, and the anomalous 1979 August 18 northeast Balapan explosion which exhibits opposite polarity, azimuth-independent source component U 1 compared to an explosion. Modelling these tests shows that Rayleigh wave amplitudes generated by the damage source are nearly as large as or larger than amplitudes from the explosion. As such, destructive interference can be drastic, introducing metastable conditions due to the sensitivity of reduced amplitudes to Rayleigh wave initial phase angles of the explosion and damage sources. This meta-stability is a likely source of scatter in M s -yield scaling observations. The agreement of observed scaling exponent b with the model suggests that the damage source strength does not vary much with yield, in contrast to explosions conducted in weak media where M s scaling rates are greater than the model predicts, and the yield dependence of the damage source strength is significant. This difference in scaling behaviour is a consequence of source medium material properties.

Description: The cancer genome is abnormal genome, and the ability to monitor its sequence had undergone a technological revolution. Yet prognosis and diagnosis remain an expert-based decision, with only limited abilities to provide machine-based decisions. We introduce a heterogeneity-based method for stratifying and visualizing whole-genome sequencing (WGS) reads. This method uses the heterogeneity within WGS reads to markedly reduce the dimensionality of next-generation sequencing data; it is available through the tool HiBS (Heterogeneity-Based Subclassification) that allows cancer sample classification. We validated HiBS using 〉200 WGS samples from nine different cancer types from The Cancer Genome Atlas (TCGA). With HiBS, we show progress with two WGS related issues: (i) differentiation between normal (NB) and tumor (TP) samples based solely on the information structure of their WGS data, and (ii) identification of specific regions of chromosomal amplification/deletion and their association with tumor stage. By comparing results to those obtained through available WGS analyses tools, we demonstrate some of the novelties obtained by the approach implemented in HiBS and also show nearly perfect normal/tumor classification, used to identify known and unknown chromosomal aberrations. Finally, the HiBS index has been associated with breast cancer tumor stage.

Description: Recent evidence suggests that many endogenous circular RNAs (circRNAs) may play roles in biological processes. However, the expression patterns and functions of circRNAs in human diseases are not well understood. Computationally identifying circRNAs from total RNA-seq data is a primary step in studying their expression pattern and biological roles. In this work, we have developed a computational pipeline named UROBORUS to detect circRNAs in total RNA-seq data. By applying UROBORUS to RNA-seq data from 46 gliomas and normal brain samples, we detected thousands of circRNAs supported by at least two read counts, followed by successful experimental validation on 24 circRNAs from the randomly selected 27 circRNAs. UROBORUS is an efficient tool that can detect circRNAs with low expression levels in total RNA-seq without RNase R treatment. The circRNAs expression profiling revealed more than 476 circular RNAs differentially expressed in control brain tissues and gliomas. Together with parental gene expression, we found that circRNA and its parental gene have diversified expression patterns in gliomas and control brain tissues. This study establishes an efficient and sensitive approach for predicting circRNAs using total RNA-seq data. The UROBORUS pipeline can be accessed freely for non-commercial purposes at http://uroborus.openbioinformatics.org/ .

Description: Recent studies show that RNA-binding proteins (RBPs) and microRNAs (miRNAs) function in coordination with each other to control post-transcriptional regulation (PTR). Despite this, the majority of research to date has focused on the regulatory effect of individual RBPs or miRNAs. Here, we mapped both RBP and miRNA binding sites on human 3'UTRs and utilized this collection to better understand PTR. We show that the transcripts that lack competition for HuR binding are destabilized more after HuR depletion. We also confirm this finding for PUM1(2) by measuring genome-wide expression changes following the knockdown of PUM1(2) in HEK293 cells. Next, to find potential cooperative interactions, we identified the pairs of factors whose sites co-localize more often than expected by random chance. Upon examining these results for PUM1(2), we found that transcripts where the sites of PUM1(2) and its interacting miRNA form a stem-loop are more stabilized upon PUM1(2) depletion. Finally, using dinucleotide frequency and counts of regulatory sites as features in a regression model, we achieved an AU-ROC of 0.86 in predicting mRNA half-life in BEAS-2B cells. Altogether, our results suggest that future studies of PTR must consider the combined effects of RBPs and miRNAs, as well as their interactions.

Description: Annotation of protein-coding genes is very important in bioinformatics and biology and has a decisive influence on many downstream analyses. Homology-based gene prediction programs allow for transferring knowledge about protein-coding genes from an annotated organism to an organism of interest. Here, we present a homology-based gene prediction program called GeMoMa. GeMoMa utilizes the conservation of intron positions within genes to predict related genes in other organisms. We assess the performance of GeMoMa and compare it with state-of-the-art competitors on plant and animal genomes using an extended best reciprocal hit approach. We find that GeMoMa often makes more precise predictions than its competitors yielding a substantially increased number of correct transcripts. Subsequently, we exemplarily validate GeMoMa predictions using Sanger sequencing. Finally, we use RNA-seq data to compare the predictions of homology-based gene prediction programs, and find again that GeMoMa performs well. Hence, we conclude that exploiting intron position conservation improves homology-based gene prediction, and we make GeMoMa freely available as command-line tool and Galaxy integration.

Description: CircRNAs are novel members of the non-coding RNA family. For several decades circRNAs have been known to exist, however only recently the widespread abundance has become appreciated. Annotation of circRNAs depends on sequencing reads spanning the backsplice junction and therefore map as non-linear reads in the genome. Several pipelines have been developed to specifically identify these non-linear reads and consequently predict the landscape of circRNAs based on deep sequencing datasets. Here, we use common RNAseq datasets to scrutinize and compare the output from five different algorithms; circRNA_finder, find_circ, CIRCexplorer, CIRI, and MapSplice and evaluate the levels of bona fide and false positive circRNAs based on RNase R resistance. By this approach, we observe surprisingly dramatic differences between the algorithms specifically regarding the highly expressed circRNAs and the circRNAs derived from proximal splice sites. Collectively, this study emphasizes that circRNA annotation should be handled with care and that several algorithms should ideally be combined to achieve reliable predictions.

Description: We present a systematic study of seismic velocity changes associated with a megathrust and five strong crustal earthquakes in Japan. We perform both cross-correlation and single-station cross-correlation analysis for station pairs and stations, respectively. The correlation of ambient seismic noise allows us to reconstruct the Green's functions of the wave propagation. By relating the coda parts of the daily Green's functions with the long-term reference Green's functions, shear wave velocity changes are determined. We analyse data from four areas in Japan where large earthquakes occurred: Iwate-Miyagi (2008 M W 6.9 Iwate-Miyagi Nairiku earthquake), Niigata (2004 M W 6.6 Chūetsu, 2007 M W 6.6 Chūetsu-oki and 2011 M W 6.2 Nagano/Niigata earthquakes), Noto Peninsula (2007 M W 6.7 Noto Hantō earthquake) and Fukuoka (2005 M W 6.6 Fukuoka earthquake). In all areas, we analyse time-series which start before the respective earthquakes and last until after the 2011 M W 9.0 Tōhoku-oki earthquake. The analysis in five different frequency ranges between 0.125 and 4.0 Hz yields time-series of the velocity changes for the different station pairs or stations. At the time of the respective earthquakes, we observe coseismic velocity drops in all areas which are followed by a partial post-seismic recovery process. For the Tōhoku-oki earthquake, coseismic velocity drops can also be observed in all regions. There is a general trend of increasing coseismic velocity drops with frequency in all four areas. The largest coseismic drops are observed close to the fault zones. Over the observed time range, the post-seismic recovery is only partial and around half of the coseismic velocity drops do not recover. The characteristic recovery times for the recovering part are similar in all areas and frequency ranges, with an average value of 0.55 yr. We model the volumetric strain changes for the different earthquakes and find that the observed pattern of the coseismic velocity drops cannot be explained by these models. The coseismic velocity drops at the different stations are better related with the peak ground velocities and the associated dynamic strain than with the peak ground accelerations, but the correlation is still poor. This suggests that non-linear effects caused by the strong ground motion during the earthquake can explain at least part of the coseismic velocity drops.

Description: Subducting slabs create strong short wavelength seismic anomalies in the upper mantle where much of Earth's seismicity is located. As such, they have the potential to bias longer wavelength seismic tomography models. To evaluate the effect of subducting slabs in global tomography, we performed a series of inversions using a global synthetic shear wave traveltime data set for a theoretical slab model based on predicted thermal anomalies within slabs. The spectral element method was applied to predict the traveltime anomalies produced by the 3-D slab model for paths corresponding to our current data used in actual tomography models. Inversion tests have been conducted first using the raw traveltime anomalies to check how well the slabs can be imaged in global tomography without the effect of earthquake mislocation. Our results indicate that most of the slabs can be identified in the inversion result but with smoothed and reduced amplitude. The recovery of the total mass anomaly in slab regions is about 88 per cent. We then performed another inversion test to investigate the effect of mislocation caused by subducting slabs. We found that source mislocation largely removes slab signal and significantly degrades the imaging of subducting slabs—potentially reducing the recovery of mass anomalies in slab regions to only 41 per cent. We tested two source relocation procedures—an iterative relocation inversion and joint relocation inversion. Both methods partially recover the true source locations and improve the inversion results, but the joint inversion method worked significantly better than the iterative method. In all of our inversion tests, the amplitudes of artefact structures in the lower mantle caused by the incorrect imaging of slabs (up to ~0.5 per cent S velocity anomalies) are comparable to some large-scale lower-mantle heterogeneities seen in global tomography studies. Based on our inversion tests, we suggest including a-priori subducting slabs in the starting models in global tomography studies and use joint relocation in the inversion.

Description: The crust at the junction of Qinling Orogen, Yangtze Craton and NE Tibetan Plateau bears imprints of the Triassic collision and later intracontinental orogeny between the Qinling Orogen and the Yangtze Craton, and the Cenozoic growth of Tibetan Plateau. Investigating detailed crustal structures in this region helps to better understand these tectonic processes. In this study, we construct a 3-D crustal Vs model using seismic ambient noise data recorded at 321 seismic stations. Ambient noise tomography is performed to generate Rayleigh wave phase velocity maps at 8–50 s periods, which are then inverted for a 3-D isotropic Vs model using a Bayesian Monte Carlo method. Our 3-D model reveals deep-rooted high velocities beneath the Hannan-Micang and Shennong-Huangling Domes, which are located on the west and east sides of the Dabashan Orocline. Similar high velocities are observed in the upper/mid crust of the western Qinling Orogen. We suggest the crustal-scale bodies with high velocity beneath the two domes and the western Qinling Orogen may represent mechanically strong rocks, which not only assisted the formation of the major Dabashan Orocline during late Mesozoic intracontinental orogeny, but also have impeded the northeastward expansion of the Tibetan Plateau during the Cenozoic era.

Description: The collective excitation of city structures by a seismic wavefield and the subsequent multiple Structure-Soil-Structure Interactions (SSSIs) between the buildings are usually disregarded in conventional seismology and earthquake engineering practice. The objective here is to qualify and quantify these complex multiple SSSIs through the design of an elementary study case, which serves as a benchmark for theoretical, numerical and experimental crossed-analysis. The experimental specimen consists of an idealized site-city setup with up to 37 anisotropic resonant structures arranged at the top surface of an elastic layer and in co-resonance with it. The experimental data from shaking table measurements is compared with the theoretical and numerical results provided respectively by an equivalent city-impedance model derived analytically from homogenization in the long-wavelength approximation and a model based on boundary elements. The signatures of the site-city interactions are identified in the frequency, time and space domain, and in particular consist of a frequency-dependent free/rigid switch in the surface condition at the city resonance, beatings in the records and the depolarization of the wavefield. A parametric study on the city density shows that multiple SSSIs among the city structures (five are sufficient) can have significant effects on both the seismic response of its implantation site and that of the buildings. Key parameters are provided to assess site-city interactions in the low seismic frequency range: They involve the mass and rigidity of the city compared to those of the soil and the damping of the building.

Description: We present a novel 3-D pre-stack Kirchhoff depth migration (PKDM) method for teleseismic receiver functions. The proposed algorithm considers the effects of diffraction, scattering and traveltime alteration caused by 3-D volumetric heterogeneities. It is therefore particularly useful for imaging complex 3-D structures such as dipping discontinuities, which is hard to accomplish with traditional methods. The scheme is based on the acoustic wave migration principle, where at each time step of the receiver function, the energy is migrated back to the ensemble of potential conversion points in the image, given a smooth 3-D reference model. Traveltimes for P and S waves are computed with an efficient eikonal solver, the fast marching method. We also consider elastic scattering patterns, where the amplitude of converted S waves depends on the angle between the incident P wave and the scattered S wave. Synthetic experiments demonstrate the validity of the method for a variety of dipping angle discontinuities. Comparison with the widely used common conversion point (CCP) stacking method reveals that our migration shows considerable improvement. For example, the effect of multiple reflections that usually produce apparent discontinuities is avoided. The proposed approach is practical, computationally efficient, and is therefore a potentially powerful alternative to standard CCP methods for imaging large-scale continental structure under dense networks.

Description: Wavelet estimation is an essential step in qualitatively and quantitatively analysing and interpreting seismic data. Applications span from seismic data quality assessment to well ties and seismic inversion. Wavelet estimation methods can be roughly separated into two approaches, data driven inversion methods and analytical definitions. We present a new analytical wavelet definition, which is based on Hermite basis functions. This wavelet model contains four parameters, which correspond to wavelet magnitude, phase, wavelet length and bandwidth. One of our main motivations for this development was to define a compact wavelet representation and an intrinsic parameter uncertainty assessment workflow, which allows us to quantify uncertainties in estimated wavelets, as well as the generation of wavelet realizations to be used, for example, in statistical seismic amplitude inversions. We present a statistical workflow to estimate the model parameters and to explore their posterior uncertainties given well log data and seismic amplitude data. This includes sampling the posterior distribution of the four wavelet parameters using Markov Chain Monte Carlo methods. We then discuss the applicability, limitations and challenges of the approach with the help of synthetic data and a North Sea data set with well logs and processed seismic amplitudes, where we also compare our method to Bayesian least-squares and a commercial wavelet estimation routine. Realizations of wavelets based on the optimized parameters and their uncertainties appear to sample the wavelet space well with reasonable variations in wavelet length, phase and amplitude while not introducing random fluctuations or wavelet lobes. The results indicate that the compact wavelet representation allows for an efficient and rather stable wavelet estimation workflow that achieves useful results in the presence of noisy data.

Description: We investigate the upper crustal velocity structure beneath the Tehri region of the Garhwal Himalaya. The investigated region is situated within the 700-km-long central seismic gap of the Himalaya that has experienced three gap-filling earthquakes since 1991 including the recent 2015 Nepal earthquake ( M w 7.8). The local tomographic inversion is based on a data set of 1365 events collected from 2008 January to 2012 December by a 12-station local network that covers an area of about 100  x  80 km around Tehri Dam. We perform a simultaneous inversion for P- and S -wave velocity anomalies. Tomograms are interpreted in the backdrop of the regional geological and tectonic framework of the region. The spatial distribution of relocated events from the 3-D velocity model has shed new light on the pattern of seismicity in the vicinity of the Main Central thrust (MCT), and has elucidated the structure of the underthrusting Indian plate. Our model exhibits a significant negative velocity anomaly up to ~5 per cent beneath the central part of the Garhwal Inner Lesser Himalaya, and a P -wave low velocity anomaly near the Chamoli region. The seismicity zone around the Chamoli region may be attributed to the presence of fluid-filled rocks. Furthermore, an area with ~3–4 per cent positive velocity anomaly is delineated to the northwest of the Uttarkashi thrust in the vicinity of the MCT. Significant findings of the study include: a flat-ramp-flat-type subsurface geometry of the underthrusting Indian plate below the Garhwal Himalaya, high-velocity images representing the trend and configuration of Delhi–Haridwar ridge below the Sub Himalaya and Lesser Himalaya and a seismically active zone representing geometrical asperity on the basement thrust in the vicinity of the MCT.

Description: Background: Fusion transcripts are formed by either fusion genes (DNA level) or trans-splicing events (RNA level). They have been recognized as a promising tool for diagnosing, subtyping and treating cancers. RNA-seq has become a precise and efficient standard for genome-wide screening of such aberration events. Many fusion transcript detection algorithms have been developed for paired-end RNA-seq data but their performance has not been comprehensively evaluated to guide practitioners. In this paper, we evaluated 15 popular algorithms by their precision and recall trade-off, accuracy of supporting reads and computational cost. We further combine top-performing methods for improved ensemble detection. Results: Fifteen fusion transcript detection tools were compared using three synthetic data sets under different coverage, read length, insert size and background noise, and three real data sets with selected experimental validations. No single method dominantly performed the best but SOAPfuse generally performed well, followed by FusionCatcher and JAFFA. We further demonstrated the potential of a meta-caller algorithm by combining top performing methods to re-prioritize candidate fusion transcripts with high confidence that can be followed by experimental validation. Conclusion: Our result provides insightful recommendations when applying individual tool or combining top performers to identify fusion transcript candidates.

Description: Slow-slip events (SSEs) are common at subduction zone faults where large mega earthquakes occur. We report here that one of the best-recorded moderate size continental earthquake, the 2009 April 6 moment magnitude ( M w ) 6.3 L'Aquila (Italy) earthquake, was preceded by a 5.9 M w SSE that originated from the decollement beneath the reactivated normal faulting system. The SSE is identified from a rigorous analysis of continuous GPS stations and occurred on the 12 February and lasted for almost two weeks. It coincided with a burst in the foreshock activity with small repeating earthquakes migrating towards the main-shock hypocentre as well as with a change in the elastic properties of rocks in the fault region. The SSE has caused substantial stress loading at seismogenic depths where the magnitude 4.0 foreshock and M w 6.3 main shock nucleated. This stress loading is also spatially correlated with the lateral extent of the aftershock sequence.

Description: For more than two decades, the number of data and model parameters in seismic tomography problems has exceeded the available computational resources required for application of direct computational methods, leaving iterative solvers the only option. One disadvantage of the iterative techniques is that the inverse of the matrix that defines the system is not explicitly formed, and as a consequence, the model resolution and covariance matrices cannot be computed. Despite the significant effort in finding computationally affordable approximations of these matrices, challenges remain, and methods such as the checkerboard resolution tests continue to be used. Based upon recent developments in sparse algorithms and high-performance computing resources, we show that direct methods are becoming feasible for large seismic tomography problems. We demonstrate the application of QR factorization in solving the regional P -wave structure and computing the full resolution matrix with 267 520 model parameters.

Description: Inhomogeneous noise sources surrounding stations produce asymmetric amplitudes in cross-correlation functions that yield preferential source directions. Here we show that preprocessing biases the dominant source direction estimate towards the source producing long-duration signals by down-weighting high-amplitude signals. Tests with both synthetic data and observations show that conventional preprocessing, where only earthquakes and local transients (e.g. trawling, fish impacts) are removed, is more sensitive to coherent energy, while one-bit preprocessing and running-absolute-mean preprocessing are more influenced by signal duration. Comparisons between different preprocessing methods are made on data from the Cascadia Initiative ocean bottom seismometer array, where we find that the total energy arriving from pelagic and coastal areas is similar. Moreover, pelagic-generated signals tend to be weaker but have longer duration, in contrast to coastal-generated signals that tend to be stronger but have shorter duration.

Description: Understanding earthquake clustering in space and time is important but also challenging because of complexities in earthquake patterns and the large and diverse nature of earthquake catalogues. Swarms are of particular interest because they likely result from physical changes in the crust, such as slow slip or fluid flow. Both swarms and clusters resulting from aftershock sequences can span a wide range of spatial and temporal scales. Here we test and implement a new method to identify seismicity clusters of varying sizes and discriminate them from randomly occurring background seismicity. Our method searches for the closest neighbouring earthquakes in space and time and compares the number of neighbours to the background events in larger space/time windows. Applying our method to California's San Jacinto Fault Zone (SJFZ), we find a total of 89 swarm-like groups. These groups range in size from 0.14 to 7.23 km and last from 15 min to 22 d. The most striking spatial pattern is the larger fraction of swarms at the northern and southern ends of the SJFZ than its central segment, which may be related to more normal-faulting events at the two ends. In order to explore possible driving mechanisms, we study the spatial migration of events in swarms containing at least 20 events by fitting with both linear and diffusion migration models. Our results suggest that SJFZ swarms are better explained by fluid flow because their estimated linear migration velocities are far smaller than those of typical creep events while large values of best-fitting hydraulic diffusivity are found.

Description: High-frequency seismic ground motion (10+ Hz), as needed for earthquake engineering design purposes, is largely controlled by the metre-scale structure of the earth's crust. However, the state-of-the-art velocity models poorly resolve small-scale features of the subsurface velocity and density variation. We invert 35 sonic logs (up to 3000 m in depth) in and near Los Angeles basin, CA, to obtain a statistical description of the small-scale heterogeneities of the basin. Assuming a von Karman autocorrelation function, our analysis finds that Hurst numbers, , between 0.0 and 0.2, vertical correlation lengths, a z , of 15–150 m and standard deviations of about 5 per cent characterize the variability in the borehole data. We report average parameters for Los Angeles basin of = 0.064 (0.058, 0.069) ± 0.01 (0.006, 0.012) and a z = 54 (51.1, 57.6) ± 5.9 (1.79, 9.53) m with 95 per cent confidence intervals listed in the parentheses. Despite the large depth range of the logs, there is no significant variation of the statistical parameters with depth. Our analysis of 371 depth-averaged shear wave velocities in the upper 30 m, V s 30 , provides only an upper bound of basin scale-length estimates due to the coarse sampling distance, with a Hurst number of about 0.3 and lateral correlation lengths, a x , of 5–10 km.

Description: A subsurface extended image gather obtained during controlled-source depth imaging yields a blurred kernel of an interface reflection operator. This reflectivity kernel or reflection function is comprised of the interface plane-wave reflection coefficients and so, in principle, the gather contains amplitude versus offset or angle information. We present a modelling theory for extended image gathers that accounts for variable illumination and blurring, under the assumption of a good migration-velocity model. The method involves forward modelling as well as migration or back propagation so as to define a receiver-side blurring function, which contains the effects of the detector array for a given shot. Composition with the modelled incident wave and summation over shots then yields an overall blurring function that relates the reflectivity to the extended image gather obtained from field data. The spatial evolution or instability of blurring functions is a key concept and there is generally not just spatial blurring in the apparent reflectivity, but also slowness or angle blurring. Gridded blurring functions can be estimated with, for example, a reverse-time migration modelling engine. A calibration step is required to account for ad hoc band limitedness in the modelling and the method also exploits blurring-function reciprocity. To demonstrate the concepts, we show numerical examples of various quantities using the well-known SIGSBEE test model and a simple salt-body overburden model, both for 2-D. The moderately strong slowness/angle blurring in the latter model suggests that the effect on amplitude versus offset or angle analysis should be considered in more realistic structures. Although the description and examples are for 2-D, the extension to 3-D is conceptually straightforward. The computational cost of overall blurring functions implies their targeted use for the foreseeable future, for example, in reservoir characterization. The description is for scalar waves, but the extension to elasticity is foreseeable and we emphasize the separation of the overburden and survey-geometry blurring effects from the nature of the target scatterer.

Description: Regulatory DNA elements, short genomic segments that regulate gene expression, have been implicated in developmental disorders and human disease. Despite this clinical urgency, only a small fraction of the regulatory DNA repertoire has been confirmed through reporter gene assays. The overall success rate of functional validation of candidate regulatory elements is low. Moreover, the number and diversity of datasets from which putative regulatory elements can be identified is large and rapidly increasing. We generated a flexible and user-friendly tool to integrate the information from different types of genomic datasets, e.g. ATAC-seq, ChIP-seq, conservation, aiming to increase the ease and success rate of functional prediction. To this end, we developed the EMERGE program that merges all datasets that the user considers informative and uses a logistic regression framework, based on validated functional elements, to set optimal weights to these datasets. ROC curve analysis shows that a combination of datasets leads to improved prediction of tissue-specific enhancers in human, mouse and Drosophila genomes. Functional assays based on this prediction can be expected to have substantially higher success rates. The resulting integrated signal for prediction of functional elements can be plotted in a build-in genome browser or exported for further analysis.

Description: Regulation of gene expression requires both transcription factor (TFs) and epigenetic modifications, and interplays between the two types of factors have been discovered. However study of relationships between chromatin features and TF–TF co-occupancy remains limited. Here, we revealed the relationship by first illustrating distinct profile patterns of chromatin features related to different binding events, including single TF binding and TF–TF co-occupancy of 71 TFs from five human cell lines. We further implemented statistical analyses to demonstrate the relationship by accurately predicting co-occupancy genome-widely using chromatin features including DNase I hypersensitivity, 11 histone modifications (HMs) and GC content. Remarkably, our results showed that the combination of chromatin features enables accurate predictions across the five cells. For individual chromatin features, DNase I enables high and consistent predictions. H3K27ac, H3K4me 2, H3K4me3 and H3K9ac are more reliable predictors than other HMs. Although the combination of 11 HMs achieves accurate predictions, their predictive ability varies considerably when a model obtained from one cell is applied to others, indicating relationship between HMs and TF–TF co-occupancy is cell type dependent. GC content is not a reliable predictor, but the addition of GC content to any other features enhances their predictive ability. Together, our results elucidate a strong relationship between TF–TF co-occupancy and chromatin features.

Description: The spatial organization of the genome influences cellular function, notably gene regulation. Recent studies have assessed the three-dimensional (3D) co-localization of functional annotations (e.g. centromeres, long terminal repeats) using 3D genome reconstructions from Hi-C (genome-wide chromosome conformation capture) data; however, corresponding assessments for continuous functional genomic data (e.g. chromatin immunoprecipitation-sequencing (ChIP-seq) peak height) are lacking. Here, we demonstrate that applying bump hunting via the patient rule induction method (PRIM) to ChIP-seq data superposed on a Saccharomyces cerevisiae 3D genome reconstruction can discover ‘functional 3D hotspots’, regions in 3-space for which the mean ChIP-seq peak height is significantly elevated. For the transcription factor Swi6, the top hotspot by P -value contains MSB2 and ERG11 – known Swi6 target genes on different chromosomes. We verify this finding in a number of ways. First, this top hotspot is relatively stable under PRIM across parameter settings. Second, this hotspot is among the top hotspots by mean outcome identified by an alternative algorithm, k -Nearest Neighbor ( k -NN) regression. Third, the distance between MSB2 and ERG11 is smaller than expected (by resampling) in two other 3D reconstructions generated via different normalization and reconstruction algorithms. This analytic approach can discover functional 3D hotspots and potentially reveal novel regulatory interactions.

Description: We analyse the influence of the solid Earth tides and ocean loading on the occurrence time of Southern California earthquakes. For each earthquake, we calculate tidal Coulomb failure stress and stress rate on a fault plane that is assumed to be controlled by the orientation of the adjacent fault. To reduce bias when selecting data for testing the tide-earthquake relationship, we create four earthquake catalogues containing events within 1, 1.5, 2.5 and 5 km of nearest faults. We investigate the difference in seismicity rates at times of positive and negative tidal stresses/stress rates given three different cases. We consider seismicity rates during times of positive versus negative stress and stress rate, as well as 2 and 3 hr surrounding the local tidal stress extremes. We find that tidal influence on earthquake occurrence is found to be statistically non-random only in close proximity to tidal extremes meaning that magnitude of tidal stress plays an important role in tidal triggering. A non-random tidal signal is observed for the reverse events. Along with a significant increase in earthquake rates around tidal Coulomb stress maxima, the strength of tidal correlation is found to be closely related to the amplitude of the peak tidal Coulomb stress ( p ). The most effective tidal triggering is found for p ≥ 1 kPa, which is much smaller than thresholds suggested for static and dynamic triggering of aftershocks.

Description: Hidden Markov models (HMMs) have been extensively used to dissect the genome into functionally distinct regions using data such as RNA expression or DNA binding measurements. It is a challenge to disentangle processes occurring on complementary strands of the same genomic region. We present the double-stranded HMM (dsHMM), a model for the strand-specific analysis of genomic processes. We applied dsHMM to yeast using strand specific transcription data, nucleosome data, and protein binding data for a set of 11 factors associated with the regulation of transcription.The resulting annotation recovers the mRNA transcription cycle (initiation, elongation, termination) while correctly predicting strand-specificity and directionality of the transcription process. We find that pre-initiation complex formation is an essentially undirected process, giving rise to a large number of bidirectional promoters and to pervasive antisense transcription. Notably, 12% of all transcriptionally active positions showed simultaneous activity on both strands. Furthermore, dsHMM reveals that antisense transcription is specifically suppressed by Nrd1, a yeast termination factor.

Description: Chromatin immunoprecipitation with massively parallel sequencing (ChIP-seq) is widely used to identify binding sites for a target protein in the genome. An important scientific application is to identify changes in protein binding between different treatment conditions, i.e. to detect differential binding. This can reveal potential mechanisms through which changes in binding may contribute to the treatment effect. The csaw package provides a framework for the de novo detection of differentially bound genomic regions. It uses a window-based strategy to summarize read counts across the genome. It exploits existing statistical software to test for significant differences in each window. Finally, it clusters windows into regions for output and controls the false discovery rate properly over all detected regions. The csaw package can handle arbitrarily complex experimental designs involving biological replicates. It can be applied to both transcription factor and histone mark datasets, and, more generally, to any type of sequencing data measuring genomic coverage. csaw performs favorably against existing methods for de novo DB analyses on both simulated and real data. csaw is implemented as a R software package and is freely available from the open-source Bioconductor project.

Description: Incremental selection within a population, defined as limited fitness changes following mutation, is an important aspect of many evolutionary processes. Strongly advantageous or deleterious mutations are detected using the synonymous to non-synonymous mutations ratio. However, there are currently no precise methods to estimate incremental selection. We here provide for the first time such a detailed method and show its precision in multiple cases of micro-evolution. The proposed method is a novel mixed lineage tree/sequence based method to detect within population selection as defined by the effect of mutations on the average number of offspring. Specifically, we propose to measure the log of the ratio between the number of leaves in lineage trees branches following synonymous and non-synonymous mutations. The method requires a high enough number of sequences, and a large enough number of independent mutations. It assumes that all mutations are independent events. It does not require of a baseline model and is practically not affected by sampling biases. We show the method's wide applicability by testing it on multiple cases of micro-evolution. We show that it can detect genes and inter-genic regions using the selection rate and detect selection pressures in viral proteins and in the immune response to pathogens.

Description: We use seismic waves for a magnitude 6.3 intermediate-depth (160 km) earthquake in the Bucaramanga Nest, Colombia, to infer a complex rupture process with two distinct stages, characterized by different rupture velocities possibly controlled by the evolution of strength on the fault. Our integrated data processing permitted to precisely characterize the multistage rupture and the presence of a strong weakening event. The resulting seismic radiation is interpreted as resulting from an extreme weakening due to a cascading thermal shear runaway, with an initial inefficient radiation process followed by a fast and dynamic efficient rupture. Our results imply dynamic complexity of the seismic rupture deep inside the Earth, and may help to give some new insights about the physical mechanism of intermediate-depth earthquakes.

Description: We investigated the shallow structure of the Solfatara, a volcano within the Campi Flegrei caldera, southern Italy, using surface waves as a diagnostic tool. We analysed data collected during the RICEN campaign, where a 3-D active seismic experiment was performed on a dense regular grid of 90 m  x  115 m using a Vibroseis as the seismic source. After removal of the source time function, we analysed the surface wave contribution to the Green's function. Here, a 1-D approximation can hold for subgrids of 40 m  x  40 m. Moreover, we stacked all of the signals in the subgrid according to source–receiver distance bins, despite the absolute location of the source and the receiver, to reduce the small-scale variability in the data. We then analysed the resulting seismic sections in narrow frequency bands between 7 and 25 Hz. We obtained phase and group velocities from a grid search, and a cost function based on the spatial coherence of both the waveforms and their envelopes. We finally jointly inverted the dispersion curves of the phase and group velocities to retrieve a 1-D S -wave model local to the subgrid. Together, the models provided a 3-D description of the S -wave model in the area. We found that the maximum penetration depth is 15 m. In the first 4 m, we can associate the changes in the S -wave field to the temperature gradient, while at greater depths, the seismic images correlate with the resistivity maps, which indicate the water layer close to the Fangaia area and an abrupt variation moving towards the northeast.

Description: An important challenge in cancer genomics is precise detection of structural variations (SVs) by high-throughput short-read sequencing, which is hampered by the high false discovery rates of existing analysis tools. Here, we propose an accurate SV detection method named COSMOS, which compares the statistics of the mapped read pairs in tumor samples with isogenic normal control samples in a distinct asymmetric manner. COSMOS also prioritizes the candidate SVs using strand-specific read-depth information. Performance tests on modeled tumor genomes revealed that COSMOS outperformed existing methods in terms of F-measure. We also applied COSMOS to an experimental mouse cell-based model, in which SVs were induced by genome engineering and gamma-ray irradiation, followed by polymerase chain reaction-based confirmation. The precision of COSMOS was 84.5%, while the next best existing method was 70.4%. Moreover, the sensitivity of COSMOS was the highest, indicating that COSMOS has great potential for cancer genome analysis.

Description: Microseisms in the period of 2–10 s are generated in deep oceans and near coastal regions. It is common for microseisms from multiple sources to arrive at the same time at a given seismometer. It is therefore desirable to be able to measure multiple slowness vectors accurately. Popular ways to estimate the direction of arrival of ocean induced microseisms are the conventional (fk) or adaptive (Capon) beamformer. These techniques give robust estimates, but are limited in their resolution capabilities and hence do not always detect all arrivals. One of the limiting factors in determining direction of arrival with seismic arrays is the array response, which can strongly influence the estimation of weaker sources. In this work, we aim to improve the resolution for weaker sources and evaluate the performance of two deconvolution algorithms, Richardson–Lucy deconvolution and a new implementation of CLEAN-PSF. The algorithms are tested with three arrays of different aperture (ASAR, WRA and NORSAR) using 1 month of real data each and compared with the conventional approaches. We find an improvement over conventional methods from both algorithms and the best performance with CLEAN-PSF. We then extend the CLEAN-PSF framework to three components (3C) and evaluate 1 yr of data from the Pilbara Seismic Array in northwest Australia. The 3C CLEAN-PSF analysis is capable in resolving a previously undetected Sn phase.

Description: We investigated a detailed spatial distribution of principal stress axis orientations in the source area of the 2003 M6.4 Northern Miyagi Prefecture earthquake that occurred in the forearc of northeastern Japan. Aftershock hypocentres were precisely relocated by applying the double difference method to arrival time data obtained at temporary stations as well as at surrounding routine stations. We picked many P -wave polarity data from seismograms at these stations, which enabled us to obtain 312 well-determined focal mechanism solutions. Stress tensor inversions were performed by using these focal mechanism data. The results show that quite a lot of focal mechanisms are difficult to explain by the uniform stress field, especially near the large slip area of the main-shock rupture. Stress tensor inversions at the location of individual earthquakes show that 1 axes are orientated mainly to WSW–ENE in the northern part of the source area, while they are oriented to NW–SE in the southern part. This spatial pattern is roughly similar to those of the static stress change by the main shock, which suggests that the observed spatially heterogeneous stress field was formed by the static stress change. If this is the case, the deviatoric stress magnitude before the main shock was very small. Another possibility is the heterogeneous stress field observed after the main shock had existed even before the main shock, although we do not know why it was formed. Unfavourable orientation of the main shock fault with respect to this stress field suggests that the fault is not strong in this case too.

Description: Gassmann's fluid substitution theory is commonly used to predict seismic velocity change upon change in saturation, and is hence essential for 4-D seismic and AVO studies. This paper addresses the basics assumptions of the Gassmann theory, in order to see how well they are fulfilled in controlled laboratory experiments. Our focus is to investigate the sensitivity of shear modulus to fluid saturation, and the predictability of Gassmann's fluid substitution theory for P -wave modulus. Ultrasonic P - and S -wave velocities in dry and saturated (3.5 wt per cent NaCl) unconsolidated clean sands (Ottawa and Columbia) were measured in an oedometer test system (uniaxial strain conditions) over a range of 0.5–10 MPa external vertical stress. This study shows shear modulus hardening upon brine saturation, which is consistent with previous data found in the literature. Analysis of the data shows that most of the hardening of the ultrasonic shear modulus may be explained by Biot dispersion. Isotropic Gassmann's fluid substitution is found to underestimate the P -wave modulus upon fluid saturation. However, adding the Biot dispersion effect improves the prediction. More work is required to obtain good measurements of parameters influencing dispersion, such as tortuosity, which is very ambiguous and challenging to measure accurately.

Description: Modern tectonic studies often use regional moment tensors (RMTs) to interpret the seismotectonic framework of an earthquake or earthquake sequence; however, despite extensive use, little existing work addresses RMT parameter uncertainty. Here, we quantify how network geometry and faulting style affect RMT sensitivity. We examine how data-model fits change with fault plane geometry (strike and dip) for varying station configurations. We calculate the relative data fit for incrementally varying geometries about a best-fitting solution, applying our workflow to real and synthetic seismograms for both real and hypothetical station distributions and earthquakes. Initially, we conduct purely observational tests, computing RMTs from synthetic seismograms for hypothetical earthquakes and a series of well-behaved network geometries. We then incorporate real data and station distributions from the International Maule Aftershock Deployment (IMAD), which recorded aftershocks of the 2010 M W 8.8 Maule earthquake, and a set of regional stations capturing the ongoing earthquake sequence in Oklahoma and southern Kansas. We consider RMTs computed under three scenarios: (1) real seismic records selected for high data quality; (2) synthetic seismic records with noise computed for the observed source-station pairings and (3) synthetic seismic records with noise computed for all possible station-source pairings. To assess RMT sensitivity for each test, we observe the ‘fit falloff’, which portrays how relative fit changes when strike or dip varies incrementally; we then derive the ranges of acceptable strikes and dips by identifying the span of solutions with relative fits larger than 90 per cent of the best fit. For the azimuthally incomplete IMAD network, Scenario 3 best constrains fault geometry, with average ranges of 45° and 31° for strike and dip, respectively. In Oklahoma, Scenario 3 best constrains fault dip with an average range of 46°; however, strike is best constrained by Scenario 1, with a range of 26°. We draw two main conclusions from this study. (1) Station distribution impacts our ability to constrain RMTs using waveform time-series; however, in some tectonic settings, faulting style also plays a significant role and (2) increasing station density and data quantity (both the number of stations and the number of individual channels) does not necessarily improve RMT constraint. These results may be useful when organizing future seismic deployments (e.g. by concentrating stations in alignment with anticipated nodal planes), and in computing RMTs, either by guiding a more rigorous data selection process for input data or informing variable weighting among the selected data (e.g. by eliminating the transverse component when strike-slip mechanisms are expected).

Description: The polarity of the first motion of a seismic signal from an earthquake is an important constraint in earthquake source inversion. Microseismic events often have low signal-to-noise ratios, which may lead to difficulties estimating the correct first-motion polarities of the arrivals. This paper describes a probabilistic approach to polarity picking that can be both automated and combined with manual picking. This approach includes a quantitative estimate of the uncertainty of the polarity, improving calculation of the polarity probability density function for source inversion. It is sufficiently fast to be incorporated into an automatic processing workflow. When used in source inversion, the results are consistent with those from manual observations. In some cases, they produce a clearer constraint on the range of high-probability source mechanisms, and are better constrained than source mechanisms determined using a uniform probability of an incorrect polarity pick.

Description: Acoustic and gravity waves propagating in planetary atmospheres have been studied intensively as markers of specific phenomena such as tectonic events or explosions or as contributors to atmosphere dynamics. To get a better understanding of the physics behind these dynamic processes, both acoustic and gravity waves propagation should be modelled in a 3-D attenuating and windy atmosphere extending from the ground to the upper thermosphere. Thus, in order to provide an efficient numerical tool at the regional or global scale, we introduce a finite difference in the time domain (FDTD) approach that relies on the linearized compressible Navier–Stokes equations with a background flow (wind). One significant benefit of such a method is its versatility because it handles both acoustic and gravity waves in the same simulation, which enables one to observe interactions between them. Simulations can be performed for 2-D or 3-D realistic cases such as tsunamis in a full MSISE-00 atmosphere or gravity-wave generation by atmospheric explosions. We validate the computations by comparing them to analytical solutions based on dispersion relations in specific benchmark cases: an atmospheric explosion, and a ground displacement forcing.

Description: We measure ellipticity of teleseismic Rayleigh waves at 95 seismic stations in Northern Italy, for wave period between 10 and 110 s, using an automatic technique and a large volume of high-quality seismic recordings from over 500 global earthquakes that occurred in 2008–2014. Northern Italy includes a wide range of crustal structures, from the wide and deep Po Plain sedimentary basin to outcropping sedimentary and crystalline rocks in the Northern Apennines and Alps. It thus provides an excellent case for studying the influence of shallow earth structure on polarization of surface waves. The ellipticity measurements show excellent spatial correlation with geological features in the region, such as high ellipticity associated with regions of low seismic velocity in the Po Plain and low ellipticity values in faster, hard rock regions in the Alps and Apennine mountains. Moreover, the observed ellipticity values also relate to the thickness of the basement, as highlighted by observed differences beneath the Alps and the Apennines. Comparison between observations and predicted ellipticity from a reference crustal model of the region show substantial fit, particularly for T ~ 38 s data. Discrepancy for shorter wave period suggests that slight modifications of the model are needed, and that the ellipticity measurements could help to better constrain the shallow crustal structure of the region. Predictions for the Po Plain are larger than the observations by a factor of four or more and transition from retrograde to prograde Rayleigh wave motion at the surface for periods of T ~ 10–13 s is predicted for seismic stations in the plain. Analysis of corresponding real data indicates a possible detection of teleseismic prograde particle motion, but the weak teleseismic earthquake signals are mixed with ambient noise signals at the predicted, short, transition periods. Detection of the period of polarity inversion from the joint analysis of earthquake and ambient noise ellipticity measurements may provide further, stringent, constraints on the structure of sedimentary basins.

Description: Kawakatsu et al. and Kawakatsu introduced a new fifth parameter, , to describe transverse isotropy (TI). Considering that characterizes the incidence angle dependence of body wave phase velocities for TI models, its relevance for body wave seismology is obvious. Here, we derive expressions for partial derivatives (sensitivity kernels) of surface wave phase velocity and normal mode eigenfrequency for the new set of five parameters. The partial derivative for is about twice as large as that for the conventional , indicating that should be more readily resolved. While partial derivatives for S velocities are not so changed, those for P velocities are significantly modified; the sensitivity for anisotropic P velocities is greatly reduced. In contrary to the suggestion by Dziewonski & Anderson and Anderson & Dziewonski, there is not much control on the anisotropic P velocities. On the other hand, the significance of for long-period seismology has become clear.