ALBERT

Description: Although numerous aquatic and terrestrial species adorn their surface with items secured from their surroundings, termed decorating, the physical and environmental factors that drive this behavior are often unclear. One of the best-known examples of this phenomenon are the decorator crabs (Majoidea), with almost 75% of species known to decorate. Here, we examined patterns of decorating in a coral reef-associated majoid, Camposcia retusa , to identify what factors determine patterns of, and investment in, decorating. Observations of natural decoration patterns indicate this species primarily decorate with sponges, fleshy and filamentous algae, and detritus. Decorations were primarily distributed on the carapace and hind walking legs which may reflect exoskeleton morphology. However, decoration cover did not decline with size, as is observed in some other majoid species, suggesting the factors driving decoration investment remain consistent throughout growth stages. From behavioral experiments, we determined that nonuniform decorating is a result of active selection, with crabs preferentially decorating their hind legs and carapace, with only small items placed on the carapace. Decorating in C. retusa appears to function primarily as an antipredator response, with crabs decorating at higher rates when access to shelter is limited. This study is the first to quantify the decorating habits of C. retusa and suggests that behaviorally-mediated decorating has a primary antipredator function. This study also highlights the value of manipulative behavioral experiments as a tool for assessing the behavioral mechanisms that drive decorating in animals.

Description: Pathogens are ubiquitous in nature and typically entail major fitness costs in their hosts. These costs can be particularly important when individuals exhibit poor immune defenses, as it is often the case during early developmental stages. Hence, selection should favor parental strategies limiting the risks of pathogen exposure and infection in their offspring. In this study, we investigated 1) whether females of the European earwig Forficula auricularia avoid areas contaminated with spores of the entomopathogenic fungus Metarhizium brunneum prior to and at egg laying, as well as 2) whether spore presence entails an increase in females’ investment into both pre-hatching forms of care and clutch quantity and quality. Our results first show that females did not avoid contaminated areas prior to and at egg laying. However, females returned to their eggs faster in presence of living spores compared to UV-killed or no spores. They were also more likely to construct a nest when in presence of both living and UV-killed spores (but only in one studied population). Finally, we found that spore presence did not influence maternal investment into egg grooming, egg gathering and egg defense, as well as into clutch quantity and quality. Overall, our results demonstrate that earwig females do not avoid contaminated environments, but could mitigate the associated costs of pathogen exposure by adjusting their level of egg care. These findings emphasize the importance of pathogens in the evolution of pre-hatching parental care and, more generally, in the emergence and maintenance of family life in nature.

Description: In this paper, we examined both stated norms of sex preference and actual sex-biases in parental investment in a Tibetan pastoralist society. We collected detailed demographic data on infant mortality, infant feeding, the length of interbirth intervals, and a decision when giving gifts, to examine sex-biased parental investment. Our results indicate a mismatch between self-reported son preference and measures of actual parental investment that favor daughters. We interpret this female-biased parental investment as a possible response to daughters generating more economic resources. However, the stated sex preferences of both sexes reflect cultural norms that appear to have remained unchanged over a long period, which may reflect the importance of male roles in the past. Our behavioral measures of parental investment are those most likely to be in the control of women (such as breastfeeding and interbirth interval), so this mismatch between stated and actual investment may be especially true of women.

Description: Behavioral syndromes are widely recognized as important for ecology and evolution, but most predictions about ecological impacts are based on conceptual models and are therefore imprecise. Borrowing insights from the theory of demographic heterogeneity, we derived insights about the population-dynamic effects of behavioral syndromes. If some individuals are consistently more aggressive than others, not just in interspecific contests, but also in foraging, mating, and antipredator behavior, then population dynamics could be affected by the resulting heterogeneity in demographic rates. We modeled a population with a boldness–aggressiveness syndrome (with the individual’s trait constant through life), showing that the mortality cost of boldness causes aggressive individuals to die earlier, on average, than their nonaggressive siblings. The equilibrium frequency of the aggressive type is strongly affected by the mortality cost of boldness, but not directly by the reproductive benefit of aggressiveness. Introducing aggressive types into a homogeneous nonaggressive population increases the average per-capita mortality rate at equilibrium; under many conditions, this reduces the equilibrium density. One such condition is that the reproductive benefit of aggression is frequency dependent and the population has evolved to equalize the expected fitness of the two types. Finally, if the intensity of aggressiveness can evolve, then the population is likely to evolve to an evolutionarily stable trait value under biologically reasonable assumptions. This analysis shows how a formal model can predict both how a syndrome affects population dynamics and how the population processes constrain evolution of the trait.

Description: “Motion dazzle” patterns are a form of defensive coloration suggested to prevent successful capture during motion by causing predators to misjudge the direction or speed of prey movement. Several studies have found results supporting this idea but little is known about the factors that favor the evolution of these antipredator colorations. A recent experimental study has suggested that the longitudinal striped patterns on the body of lizards can redirect attacks to the tail via the motion dazzle effect. Using a virtual predation experiment with humans and a phylogenetic comparative analysis, we show that evolution of longitudinal striped coloration is associated with prey size. Experiments showed that longitudinal stripes located at the anterior reduced lethal attacks (i.e., attacks directed to the anterior and centre ) but this benefit was greater for shorter prey. Our comparative analysis revealed a negative association between stripe occurrence and body length but no association between stripes and body width. Overall, our results suggest that the dazzle effect produced by stripes is more advantageous in shorter lizards than in longer ones and that the error induced by stripes might be distributed along the axis parallel to the prey trajectory. We discuss reasons why dazzle coloration could be associated with evolution of smaller body size in animals.

Description: The avian nest is an essential structure for offspring development. For adults, nest building entails costs in terms of time, energy, and exposure to predators and parasites. Amount and diversity of materials used for nest building depend on their availability and functionality in scenarios of sexual selection and parasitism. Green plants and feathers of different colors have been hypothesized to play key roles in offspring protection against pathogens, and we here experimentally assessed spotless starling ( Sturnus unicolor ) preferences for pigmented versus unpigmented feathers and for different green plants (aromatic vs. non-aromatic plants) as nest materials. We predicted a preferential selection of unpigmented feathers and aromatic plants according to the antimicrobial properties of these materials described in the literature. We evaluated these predictions during nest building and during egg-laying stages. As expected, starlings preferentially selected unpigmented feathers both before and during egg laying, while aromatic plants were preferentially selected only during the egg-laying stage. These results suggest that starlings prefer nest materials that enhance antimicrobial protection of their offspring. We also discuss some other, non-exclusive functions that might explain the observed preference for nest materials, especially with regard to their potential role in sexual signaling.

Description: There is ongoing debate about how and why the menstrual cycle affects women’s attraction to men. According to the dual sexuality hypothesis, women form pair-bond relationships with men who provide care but also obtain genetic benefits by biasing mating effort towards men with high-fitness genes during the fertile phase. By contrast, the commitment hypothesis proposes that attachment bonds with primary partners function to strengthen pair-bond relationships by enhancing in-pair attraction at the fertile phase, rather than extrapair attraction. We tested these hypotheses by measuring women’s daily sexual and emotional attraction towards men over the whole menstrual cycle. We employed 1) a urinary luteinizing hormone test to determine the day of ovulation, 2) a 5-part classification of menstrual cycle that identifies a distinct peri-ovulatory phase, and 3) individualized phase identification for each participant. There was a mid-cycle rise in extrapair sexual desire. Women gave and received more care from partners during the menstrual than the mid-cycle phases. Partner’s sexual attractiveness and mutual commitment did not moderate these findings. The results do not support either the dual sexuality or commitment hypotheses, and imply that female self-reported sexual desire is not strictly dependent on cyclic hormonal changes. Our results are more consistent with a recently proposed `spandrel’ hypothesis, positing cycle phase effects as a nonfunctional by-product of raised estradiol. Additionally, we found that, with the date of ovulation estimated by luteinizing hormone tests, 45% of ovulations were misclassified by the backward counting method, which urges caution in interpreting results based on counting methods.

Description: Corrigendum for “ Social conflict and costs of cooperation in meerkats are reflected in measures of stress hormones ”

Description: Kleptoparasitism (food stealing) is a significant behavior for animals that forage in social groups as it permits some individuals to obtain resources while avoiding the costs of searching for their own food. Evolutionary game theory has been used to model kleptoparasitism, with a series of differential equation-based compartmental models providing significant theoretical insights into behavior in kleptoparasitic populations. In this paper, we apply this compartmental modeling approach to kleptoparasitic behavior in a real foraging population of urban gulls ( Laridae ). Field data was collected on kleptoparasitism and a model developed that incorporated the same kleptoparasitic and defensive strategies available to the study population. Two analyses were conducted: 1) An assessment of whether the density of each behavior in the population was at an equilibrium. 2) An investigation of whether individual foragers were using Evolutionarily Stable Strategies in the correct environmental conditions. The results showed the density of different behaviors in the population could be at an equilibrium at plausible values for handling time and fight duration. Individual foragers used aggressive kleptoparasitic strategies effectively in the correct environmental conditions but some individuals in those same conditions failed to defend food items. This was attributed to the population being composed of 3 species that differed in competitive ability. These competitive differences influenced the strategies that individuals were able to use. Rather than gulls making poor behavioral decisions these results suggest a more complex 3-species model is required to describe the behavior of this population.

Description: In animal species in which parents provide food to their dependent young, offspring often display conspicuous begging signals. These solicitation behaviors are important components of parent–offspring communication, but it is currently unclear how they and the parental response covary with offspring dependency on parental food provisioning across species. Burying beetles ( Nicrophorus ) are well known for providing elaborate biparental care, including provisioning of begging larvae. By using a multispecies approach, we show that larval begging intensity, as well as the time parents spend provisioning, differ greatly between individuals of the 3 species: N. orbicollis , N. pustulatus , and N. vespilloides . Our results demonstrate that the most dependent offspring of N. orbicollis invest the most time in begging, whereas the most independent offspring of N. pustulatus invest the least amount of time in begging. Thus, we suggest that begging intensity differs due to intrinsic differences in nutritional need between the species rather than because of an arbitrary divergence in begging behavior. We further show that in all 3 species, females spend significantly more time provisioning than males, although there is considerable divergence between species in the extent to which females and males contribute to the provisioning of larvae. We discuss the potential selective factors leading to this diversification of offspring begging and parental provisioning in relation to the distinct variation in offspring dependence between the 3 species.

Description: Red, orange, and yellow carotenoid-based ornaments, which are widely used as sexual signals in many birds, fish, and reptiles, are known to exhibit multidimensional chromatic variation as a result of both the concentration and relative proportions of different constituent carotenoids with differing spectral properties. This is thought to reflect intrinsic variation in signaler quality, making it a useful basis for female choice. However, whether females are able to discriminate relevant variation in carotenoid concentration and/or composition independently of each other, and of other phenotypic or behavior traits, and if so, how this mediates their choice, is poorly understood. Here, female 3-spined sticklebacks ( Gasterosteus aculeatus ) were presented with computer-animated courting males that varied exclusively in the appearance of their carotenoid-based coloration; specifically, each male’s signal provided a perceptual match for carotenoid coloration expressed by live males with known underlying carotenoid content, thereby providing a biologically-relevant signal while precluding confounding traits influencing female choice. Females were able to discriminate between prospective mates solely on the basis of perceived variation in the allocation of carotenoids to males’ sexual signals, and exhibited a strong preference for males with coloration indicative of higher concentrations of carotenoids in their signal, rather than in response to perceived variation in the relative proportion of constituent carotenoids. This has important implications for our understanding of male signaling strategies and the information content of carotenoid-based sexual signals.

Description: Motivation The traditional view of cancer evolution states that a cancer genome accumulates a sequential ordering of mutations over a long period of time. However, in recent years it has been suggested that a cancer genome may instead undergo a one-time catastrophic event, such as chromothripsis , where a large number of mutations instead occur simultaneously . A number of potential signatures of chromothripsis have been proposed. In this work, we provide a rigorous formulation and analysis of the ‘ability to walk the derivative chromosome’ signature originally proposed by Korbel and Campbell. In particular, we show that this signature, as originally envisioned, may not always be present in a chromothripsis genome and we provide a precise quantification of under what circumstances it would be present. We also propose a variation on this signature, the H/T alternating fraction , which allows us to overcome some of the limitations of the original signature. Results We apply our measure to both simulated data and a previously analyzed real cancer dataset and find that the H/T alternating fraction may provide useful signal for distinguishing genomes having acquired mutations simultaneously from those acquired in a sequential fashion. Availability and implementation An implementation of the H/T alternating fraction is available at https://bitbucket.org/oesperlab/ht-altfrac . Contact loesper@carleton.edu Supplementary information Supplementary dataSupplementary data are available at Bioinformatics online.

Description: Motivation The identification of microRNA (miRNA) target sites is important. In the past decade, dozens of computational methods have been developed to predict miRNA target sites. Despite their existence, rarely does a method consider the well-known competition and cooperation among miRNAs when attempts to discover target sites. To fill this gap, we developed a new approach called CCmiR, which takes the cooperation and competition of multiple miRNAs into account in a statistical model to predict their target sites. Results Tested on four different datasets, CCmiR predicted miRNA target sites with a high recall and a reasonable precision, and identified known and new cooperative and competitive miRNAs supported by literature. Compared with three state-of-the-art computational methods, CCmiR had a higher recall and a higher precision. Availability and implementation CCmiR is freely available at http://hulab.ucf.edu/research/projects/miRNA/CCmiR . Contact xiaoman@mail.ucf.edu or haihu@cs.ucf.edu Supplementary information Supplementary dataSupplementary data are available at Bioinformatics online.

Description: Motivation Eukaryotic chromosomes adapt a complex and highly dynamic three-dimensional (3D) structure, which profoundly affects different cellular functions and outcomes including changes in epigenetic landscape and in gene expression. Making the scenario even more complex, cancer cells harbor chromosomal abnormalities [e.g. copy number variations (CNVs) and translocations] altering their genomes both at the sequence level and at the level of 3D organization. High-throughput chromosome conformation capture techniques (e.g. Hi-C), which are originally developed for decoding the 3D structure of the chromatin, provide a great opportunity to simultaneously identify the locations of genomic rearrangements and to investigate the 3D genome organization in cancer cells. Even though Hi-C data has been used for validating known rearrangements, computational methods that can distinguish rearrangement signals from the inherent biases of Hi-C data and from the actual 3D conformation of chromatin, and can precisely detect rearrangement locations de novo have been missing. Results In this work, we characterize how intra and inter-chromosomal Hi-C contacts are distributed for normal and rearranged chromosomes to devise a new set of algorithms (i) to identify genomic segments that correspond to CNV regions such as amplifications and deletions ( HiCnv ), (ii) to call inter-chromosomal translocations and their boundaries ( HiCtrans ) from Hi-C experiments and (iii) to simulate Hi-C data from genomes with desired rearrangements and abnormalities ( AveSim ) in order to select optimal parameters for and to benchmark the accuracy of our methods. Our results on 10 different cancer cell lines with Hi-C data show that we identify a total number of 105 amplifications and 45 deletions together with 90 translocations, whereas we identify virtually no such events for two karyotypically normal cell lines. Our CNV predictions correlate very well with whole genome sequencing data among chromosomes with CNV events for a breast cancer cell line ( r  = 0.89) and capture most of the CNVs we simulate using Avesim. For HiCtrans predictions, we report evidence from the literature for 30 out of 90 translocations for eight of our cancer cell lines. Furthermore, we show that our tools identify and correctly classify relatively understudied rearrangements such as double minutes and homogeneously staining regions. Considering the inherent limitations of existing techniques for karyotyping (i.e. missing balanced rearrangements and those near repetitive regions), the accurate identification of CNVs and translocations in a cost-effective and high-throughput setting is still a challenge. Our results show that the set of tools we develop effectively utilize moderately sequenced Hi-C libraries (100–300 million reads) to identify known and de novo chromosomal rearrangements/abnormalities in well-established cancer cell lines. With the decrease in required number of cells and the increase in attainable resolution, we believe that our framework will pave the way towards comprehensive mapping of genomic rearrangements in primary cells from cancer patients using Hi-C. Availability and implementation CNV calling: https://github.com/ay-lab/HiCnv , Translocation calling: https://github.com/ay-lab/HiCtrans and Hi-C simulation: https://github.com/ay-lab/AveSim . Contact ferhatay@lji.org Supplementary information Supplementary dataSupplementary data are available at Bioinformatics online.

Description: Summary Cancer hallmarks, a concept that seeks to explain the complexity of cancer initiation and development, provide a new perspective of studying cancer signaling which could lead to a greater understanding of this complex disease. However, to the best of our knowledge, there is currently a lack of tools that support such hallmark-based study of the cancer signaling network, thereby impeding the gain of knowledge in this area. We present TROVE, an user-friendly software that facilitates hallmark annotation, visualization and analysis in cancer signaling networks. In particular, TROVE facilitates hallmark analysis specific to particular cancer types. Availability and implementation Available under the Eclipse Public License from: https://sites.google.com/site/cosbyntu/softwares/trove and https://github.com/trove2017/Trove . Contact hechua@ntu.edu.sg or assourav@ntu.edu.sg

Description: Summary High-throughput screening of the host transcriptional response to various viral infections provides a wealth of data, but utilization of microarray and next generation sequencing (NGS) data for analysis can be difficult. The Ho st T ranscriptional R esponse D ata B ase (HoTResDB), allows visitors to access already processed microarray and NGS data from non-human primate models of viral hemorrhagic fever to better understand the host transcriptional response. Availability HoTResDB is freely available at http://hotresdb.bu.edu Contact jhconnor@bu.edu

Description: Motivation Structural variation, including large deletions, duplications, inversions, translocations and other rearrangements, is common in human and cancer genomes. A number of methods have been developed to identify structural variants from Illumina short-read sequencing data. However, reliable identification of structural variants remains challenging because many variants have breakpoints in repetitive regions of the genome and thus are difficult to identify with short reads. The recently developed linked-read sequencing technology from 10X Genomics combines a novel barcoding strategy with Illumina sequencing. This technology labels all reads that originate from a small number (∼5 to 10) DNA molecules ∼50 Kbp in length with the same molecular barcode. These barcoded reads contain long-range sequence information that is advantageous for identification of structural variants. Results We present Novel Adjacency Identification with Barcoded Reads (NAIBR), an algorithm to identify structural variants in linked-read sequencing data. NAIBR predicts novel adjacencies in an individual genome resulting from structural variants using a probabilistic model that combines multiple signals in barcoded reads. We show that NAIBR outperforms several existing methods for structural variant identification—including two recent methods that also analyze linked-reads—on simulated sequencing data and 10X whole-genome sequencing data from the NA12878 human genome and the HCC1954 breast cancer cell line. Several of the novel somatic structural variants identified in HCC1954 overlap known cancer genes. Availability and implementation Software is available at compbio.cs.brown.edu/software . Contact braphael@princeton.edu Supplementary information Supplementary dataSupplementary data are available at Bioinformatics online.

Description: Motivation Cancers arise as the result of somatically acquired changes in the DNA of cancer cells. However, in addition to the mutations that confer a growth advantage, cancer genomes accumulate a large number of somatic mutations resulting from normal DNA damage and repair processes as well as carcinogenic exposures or cancer related aberrations of DNA maintenance machinery. These mutagenic processes often produce characteristic mutational patterns called mutational signatures. The decomposition of a cancer genome’s mutation catalog into mutations consistent with such signatures can provide valuable information about cancer etiology. However, the results from different decomposition methods are not always consistent. Hence, one needs to be able to not only decompose a patient’s mutational profile into signatures but also establish the accuracy of such decomposition. Results We proposed two complementary ways of measuring confidence and stability of decomposition results and applied them to analyze mutational signatures in breast cancer genomes. We identified both very stable and highly unstable signatures, as well as signatures that previously have not been associated with breast cancer. We also provided additional support for the novel signatures. Our results emphasize the importance of assessing the confidence and stability of inferred signature contributions. Availability and implementation All tools developed in this paper have been implemented in an R package, called SignatureEstimation, which is available from https://www.ncbi.nlm.nih.gov/CBBresearch/Przytycka/index.cgi\#signatureestimation . Contact wojtowda@ncbi.nlm.nih.gov or przytyck@ncbi.nlm.nih.gov Supplementary information Supplementary dataSupplementary data are available at Bioinformatics online.

Description: Motivation Brain imaging genetics, which studies the linkage between genetic variations and structural or functional measures of the human brain, has become increasingly important in recent years. Discovering the bi-multivariate relationship between genetic markers such as single-nucleotide polymorphisms (SNPs) and neuroimaging quantitative traits (QTs) is one major task in imaging genetics. Sparse Canonical Correlation Analysis (SCCA) has been a popular technique in this area for its powerful capability in identifying bi-multivariate relationships coupled with feature selection. The existing SCCA methods impose either the ℓ 1 -norm or its variants to induce sparsity. The ℓ 0 -norm penalty is a perfect sparsity-inducing tool which, however, is an NP-hard problem. Results In this paper, we propose the truncated ℓ 1 -norm penalized SCCA to improve the performance and effectiveness of the ℓ 1 -norm based SCCA methods. Besides, we propose an efficient optimization algorithms to solve this novel SCCA problem. The proposed method is an adaptive shrinkage method via tuning τ . It can avoid the time intensive parameter tuning if given a reasonable small τ . Furthermore, we extend it to the truncated group-lasso (TGL), and propose TGL-SCCA model to improve the group-lasso-based SCCA methods. The experimental results, compared with four benchmark methods, show that our SCCA methods identify better or similar correlation coefficients, and better canonical loading profiles than the competing methods. This demonstrates the effectiveness and efficiency of our methods in discovering interesting imaging genetic associations. Availability and implementation The Matlab code and sample data are freely available at http://www.iu.edu/∼shenlab/tools/tlpscca/ . Contact dulei@nwpu.edu.cn or shenli@iu.edu Supplementary information Supplementary dataSupplementary data are available at Bioinformatics online.

Description: Motivation Modelling with multiple servers that use different algorithms for docking results in more reliable predictions of interaction sites. However, the scoring and comparison of all models by an expert is time-consuming and is not feasible for large volumes of data generated by such modelling. Results Quality ASsessment of DOcking Models (QASDOM) Server is a simple and efficient tool for real-time simultaneous analysis, scoring and ranking of data sets of receptor–ligand complexes built by a range of docking techniques. This meta-server is designed to analyse large data sets of docking models and rank them by scoring criteria developed in this study. It produces two types of output showing the likelihood of specific residues and clusters of residues to be involved in receptor–ligand interactions and the ranking of models. The server also allows visualizing residues that form interaction sites in the receptor and ligand sequence and displays 3D model structures of the receptor–ligand complexes. Availability http://qasdom.eimb.ru . Contact alexei.adzhubei@eimb.ru. Supplementary information Supplementary dataSupplementary data are available at Bioinformatics online.

Description: Motivation Protein–protein interactions are vital for protein function with the average protein having between three and ten interacting partners. Knowledge of precise protein–protein interfaces comes from crystal structures deposited in the Protein Data Bank (PDB), but only 50% of structures in the PDB are complexes. There is therefore a need to predict protein–protein interfaces in silico and various methods for this purpose. Here we explore the use of a predictor based on structural features and which exploits random forest machine learning, comparing its performance with a number of popular established methods. Results On an independent test set of obligate and transient complexes, our IntPred predictor performs well (MCC = 0.370, ACC = 0.811, SPEC = 0.916, SENS = 0.411) and compares favourably with other methods. Overall, IntPred ranks second of six methods tested with SPPIDER having slightly better overall performance (MCC = 0.410, ACC = 0.759, SPEC = 0.783, SENS = 0.676), but considerably worse specificity than IntPred. As with SPPIDER, using an independent test set of obligate complexes enhanced performance (MCC = 0.381) while performance is somewhat reduced on a dataset of transient complexes (MCC = 0.303). The trade-off between sensitivity and specificity compared with SPPIDER suggests that the choice of the appropriate tool is application-dependent. Availability and implementation IntPred is implemented in Perl and may be downloaded for local use or run via a web server at www.bioinf.org.uk/intpred/ . Contact andrew@bioinf.org.uk or andrew.martin@ucl.ac.uk Supplementary information Supplementary dataSupplementary data are available at Bioinformatics online.

Description: Motivation Protein function is directly related to amino acid residue composition and the dynamics of these residues. Centrality analyses based on residue interaction networks permit to identify key residues in a protein that are important for its fold or function. Such central residues and their environment constitute suitable targets for mutagenesis experiments. Predicted flexibility and changes in flexibility upon mutation provide valuable additional information for the design of such experiments. Results We combined centrality analyses with DynaMine flexibility predictions in a Cytoscape app called RINspector. The app performs centrality analyses and directly visualizes the results on a graph of predicted residue flexibility. In addition, the effect of mutations on local flexibility can be calculated. Availability and implementation The app is publicly available in the Cytoscape app store. Contact guillaume.brysbaert@univ-lille1.fr Supplementary information Supplementary dataSupplementary data are available at Bioinformatics online.

Description: Summary Fully exploiting the wealth of data in current bacterial population genomics datasets requires synthesizing and integrating different types of analysis across millions of base pairs in hundreds or thousands of isolates. Current approaches often use static representations of phylogenetic, epidemiological, statistical and evolutionary analysis results that are difficult to relate to one another. Phandango is an interactive application running in a web browser allowing fast exploration of large-scale population genomics datasets combining the output from multiple genomic analysis methods in an intuitive and interactive manner. Availability and implementation Phandango is a web application freely available for use at www.phandango.net and includes a diverse collection of datasets as examples. Source code together with a detailed wiki page is available on GitHub at https://github.com/jameshadfield/phandango . Contact jh22@sanger.ac.uk or sh16@sanger.ac.uk

Description: Summary MetExploreViz is an open source web component that can be easily embedded in any web site. It provides features dedicated to the visualization of metabolic networks and pathways and thus offers a flexible solution to analyse omics data in a biochemical context. Availability and implementation Documentation and link to GIT code repository (GPL 3.0 license) are available at this URL: http://metexplore.toulouse.inra.fr/metexploreViz/doc/ Contact contact-metexplore@inra.fr

Description: Motivation Single cell transcriptional profiling opens up a new avenue in studying the functional role of cell-to-cell variability in physiological processes. The analysis of single cell expression profiles creates new challenges due to the distributive nature of the data and the stochastic dynamics of gene transcription process. The reconstruction of gene regulatory networks (GRNs) using single cell transcriptional profiles is particularly challenging, especially when directed gene-gene relationships are desired. Results We developed SINCERITIES (SINgle CEll Regularized Inference using TIme-stamped Expression profileS) for the inference of GRNs from single cell transcriptional profiles. We focused on time-stamped cross-sectional expression data, commonly generated from transcriptional profiling of single cells collected at multiple time points after cell stimulation. SINCERITIES recovers directed regulatory relationships among genes by employing regularized linear regression (ridge regression), using temporal changes in the distributions of gene expressions. Meanwhile, the modes of the gene regulations (activation and repression) come from partial correlation analyses between pairs of genes. We demonstrated the efficacy of SINCERITIES in inferring GRNs using in silico time-stamped single cell expression data and single cell transcriptional profiles of THP-1 monocytic human leukemia cells. The case studies showed that SINCERITIES could provide accurate GRN predictions, significantly better than other GRN inference algorithms such as TSNI, GENIE3 and JUMP3. Moreover, SINCERITIES has a low computational complexity and is amenable to problems of extremely large dimensionality. Finally, an application of SINCERITIES to single cell expression data of T2EC chicken erythrocytes pointed to BATF as a candidate novel regulator of erythroid development. Availability and implementation MATLAB and R version of SINCERITIES are freely available from the following websites: http://www.cabsel.ethz.ch/tools/sincerities.html and https://github.com/CABSEL/SINCERITIES . The single cell THP-1 and T2EC transcriptional profiles are available from the original publications ( Kouno et al. , 2013 ; Richard et al. , 2016 ). The in silico single cell data are available on SINCERITIES websites. Contact rudi.gunawan@chem.ethz.ch Supplementary information Supplementary dataSupplementary data are available at Bioinformatics online.

Description: Summary Addressing deleterious effects of noncoding mutations is an essential step towards the identification of disease-causal mutations of gene regulatory elements. Several methods for quantifying the deleteriousness of noncoding mutations using artificial intelligence, deep learning and other approaches have been recently proposed. Although the majority of the proposed methods have demonstrated excellent accuracy on different test sets, there is rarely a consensus. In addition, advanced statistical and artificial learning approaches used by these methods make it difficult porting these methods outside of the labs that have developed them. To address these challenges and to transform the methodological advances in predicting deleterious noncoding mutations into a practical resource available for the broader functional genomics and population genetics communities, we developed SNPDelScore, which uses a panel of proposed methods for quantifying deleterious effects of noncoding mutations to precompute and compare the deleteriousness scores of all common SNPs in the human genome in 44 cell lines. The panel of deleteriousness scores of a SNP computed using different methods is supplemented by functional information from the GWAS Catalog, libraries of transcription factor-binding sites, and genic characteristics of mutations. SNPDelScore comes with a genome browser capable of displaying and comparing large sets of SNPs in a genomic locus and rapidly identifying consensus SNPs with the highest deleteriousness scores making those prime candidates for phenotype-causal polymorphisms. Availability and implementation https://www.ncbi.nlm.nih.gov/research/snpdelscore/ Contact ovcharen@nih.gov Supplementary information Supplementary dataSupplementary data are available at Bioinformatics online.

Description: Motivation The selection of a single nucleotide polymorphism (SNP) using bibliographic methods can be a very time-consuming task. Moreover, a SNP selected in this way may not be easily visualized in its genomic context by a standard user hoping to correlate it with other valuable information. Here we propose a web form built on top of Circos that can assist SNP-centered screening, based on their location in the genome and the regulatory modules they can disrupt. Its use may allow researchers to prioritize SNPs in genotyping and disease studies. Results SiNoPsis is bundled as a web portal. It focuses on the different structures involved in the genomic expression of a gene, especially those found in the core promoter upstream region. These structures include transcription factor binding sites (for promoter and enhancer signals), histones and promoter flanking regions. Additionally, the tool provides eQTL and linkage disequilibrium (LD) properties for a given SNP query, yielding further clues about other indirectly associated SNPs. Possible disruptions of the aforementioned structures affecting gene transcription are reported using multiple resource databases. SiNoPsis has a simple user-friendly interface, which allows single queries by gene symbol, genomic coordinates, Ensembl gene identifiers, RefSeq transcript identifiers and SNPs. It is the only portal providing useful SNP selection based on regulatory modules and LD with functional variants in both textual and graphic modes (by properly defining the arguments and parameters needed to run Circos). Availability and implementation SiNoPsis is freely available at https://compgen.bio.ub.edu/SiNoPsis/ Contact danielboloc@gmail.com Supplementary information Supplementary dataSupplementary data are available at Bioinformatics online.

Description: Motivation Inter-residue contacts in proteins have been widely acknowledged to be valuable for protein 3 D structure prediction. Accurate prediction of long-range transmembrane inter-helix residue contacts can significantly improve the quality of simulated membrane protein models. Results In this paper, we present an updated MemBrain predictor, which aims to predict transmembrane protein residue contacts. Our new model benefits from an efficient learning algorithm that can mine latent structural features, which exist in original feature space. The new MemBrain is a two-stage inter-helix contact predictor. The first stage takes sequence-based features as inputs and outputs coarse contact probabilities for each residue pair, which will be further fed into convolutional neural network together with predictions from three direct-coupling analysis approaches in the second stage. Experimental results on the training dataset show that our method achieves an average accuracy of 81.6% for the top L /5 predictions using a strict sequence-based jackknife cross-validation. Evaluated on the test dataset, MemBrain can achieve 79.4% prediction accuracy. Moreover, for the top L /5 predicted long-range loop contacts, the prediction performance can reach an accuracy of 56.4%. These results demonstrate that the new MemBrain is promising for transmembrane protein’s contact map prediction. Availability and implementation http://www.csbio.sjtu.edu.cn/bioinf/MemBrain/ Contact hbshen@sjtu.edu.cn Supplementary information Supplementary dataSupplementary data are available at Bioinformatics online.

Description: Motivation We recently published MS2LDA, a method for the decomposition of sets of molecular fragment data derived from large metabolomics experiments. To make the method more widely available to the community, here we present ms2lda.org, a web application that allows users to upload their data, run MS2LDA analyses and explore the results through interactive visualizations. Results Ms2lda.org takes tandem mass spectrometry data in many standard formats and allows the user to infer the sets of fragment and neutral loss features that co-occur together (Mass2Motifs). As an alternative workflow, the user can also decompose a data set onto predefined Mass2Motifs. This is accomplished through the web interface or programmatically from our web service. Availability and implementation The website can be found at http://ms2lda.org , while the source code is available at https://github.com/sdrogers/ms2ldaviz under the MIT license. Contact simon.rogers@glasgow.ac.uk Supplementary information Supplementary dataSupplementary data are available at Bioinformatics online.

Description: Motivation Linkage and quantitative trait loci (QTL) maps are critical tools for the study of the genetic basis of complex traits. With the advances in sequencing technology over the past decade, linkage map densities have been increasing dramatically, while the visualization tools have not kept pace. LinkageMapView is a free add-on package written in R that produces high resolution, publication-ready visualizations of linkage and QTL maps. While there is software available to generate linkage map graphics, none are freely available, produce publication quality figures, are open source and can run on all platforms. LinkageMapView can be integrated into map building pipelines as it seamlessly incorporates output from R/qtl and also accepts simple text or comma delimited files. There are numerous options within the package to build highly customizable maps, allow for linkage group comparisons, and annotate QTL regions. Availability and implementation https://cran.r-project.org/web/packages/LinkageMapView/ Contact louellet@uncc.edu

Description: Motivation Rapid and low cost sequencing of genomes enabled widespread use of genomic data in research studies and personalized customer applications, where genomic data is shared in public databases. Although the identities of the participants are anonymized in these databases, sensitive information about individuals can still be inferred. One such information is kinship. Results We define two routes kinship privacy can leak and propose a technique to protect kinship privacy against these risks while maximizing the utility of shared data. The method involves systematic identification of minimal portions of genomic data to mask as new participants are added to the database. Choosing the proper positions to hide is cast as an optimization problem in which the number of positions to mask is minimized subject to privacy constraints that ensure the familial relationships are not revealed. We evaluate the proposed technique on real genomic data. Results indicate that concurrent sharing of data pertaining to a parent and an offspring results in high risks of kinship privacy, whereas the sharing data from further relatives together is often safer. We also show arrival order of family members have a high impact on the level of privacy risks and on the utility of sharing data. Availability and implementation https://github.com/tastanlab/Kinship-Privacy Contact erman@cs.bilkent.edu.tr or oznur.tastan@cs.bilkent.edu.tr Supplementary information Supplementary dataSupplementary data are available at Bioinformatics online.

Description: Summary Precision oncology is an approach that accounts for individual differences to guide cancer management. Omics signatures have been shown to predict clinical traits for cancer patients. However, the vast amount of omics information poses an informatics challenge in systematically identifying patterns associated with health outcomes, and no general purpose data mining tool exists for physicians, medical researchers and citizen scientists without significant training in programming and bioinformatics. To bridge this gap, we built the Omics AnalySIs System for PRecision Oncology (OASISPRO), a web-based system to mine the quantitative omics information from The Cancer Genome Atlas (TCGA). This system effectively visualizes patients’ clinical profiles, executes machine-learning algorithms of choice on the omics data and evaluates the prediction performance using held-out test sets. With this tool, we successfully identified genes strongly associated with tumor stage, and accurately predicted patients’ survival outcomes in many cancer types, including adrenocortical carcinoma. By identifying the links between omics and clinical phenotypes, this system will facilitate omics studies on precision cancer medicine and contribute to establishing personalized cancer treatment plans. Availability and implementation This web-based tool is available at http://tinyurl.com/oasispro ; source codes are available at http://tinyurl.com/oasisproSourceCode . Contact khyu@stanford.edu or mpsnyder@stanford.edu Supplementary information Supplementary dataSupplementary data are available at Bioinformatics online.

Description: Summary Analysis of RNA sequencing (RNA-Seq) data have highlighted the fact that most genes undergo alternative splicing (AS) and that these patterns are tightly regulated. Many of these events are complex, resulting in numerous possible isoforms that quickly become difficult to visualize, interpret and experimentally validate. To address these challenges we developed MAJIQ-SPEL, a web-tool that takes as input local splicing variations (LSVs) quantified from RNA-Seq data and provides users with visualization and quantification of gene isoforms associated with those. Importantly, MAJIQ-SPEL is able to handle both classical (binary) and complex, non-binary, splicing variations. Using a matching primer design algorithm it also suggests to users possible primers for experimental validation by RT-PCR and displays those, along with the matching protein domains affected by the LSV, on UCSC Genome Browser for further downstream analysis. Availability and implementation Program and code will be available at http://majiq.biociphers.org/majiq-spel . Contact yosephb@upenn.edu Supplementary information Supplementary dataSupplementary data are available at Bioinformatics online.

Description: Motivation Accurate molecular structure of the protein dimer representing the elementary building block of intermediate filaments (IFs) is essential towards the understanding of the filament assembly, rationalizing their mechanical properties and explaining the effect of disease-related IF mutations. The dimer contains a ∼300-residue long α-helical coiled coil which cannot be assessed by either direct experimental structure determination or modelling using standard approaches. At the same time, coiled coils are well-represented in structural databases. Results Here we present CCFold, a generally applicable threading-based algorithm which produces coiled-coil models from protein sequence only. The algorithm is based on a statistical analysis of experimentally determined structures and can handle any hydrophobic repeat patterns in addition to the most common heptads. We demonstrate that CCFold outperforms general-purpose computational folding in terms of accuracy, while being faster by orders of magnitude. By combining the CCFold algorithm and Rosetta folding we generate representative dimer models for all IF protein classes. Availability and implementation The source code is freely available at https://github.com/biocryst/IF ; a web server to run the program is at http://pharm.kuleuven.be/Biocrystallography/cc . Contact sergei.strelkov@kuleuven.be Supplementary information Supplementary dataSupplementary data are available at Bioinformatics online.

Description: Motivation Cells process information, in part, through transcription factor (TF) networks, which control the rates at which individual genes produce their products. A TF network map is a graph that indicates which TFs bind and directly regulate each gene. Previous work has described network mapping algorithms that rely exclusively on gene expression data and ‘integrative’ algorithms that exploit a wide range of data sources including chromatin immunoprecipitation sequencing (ChIP-seq) of many TFs, genome-wide chromatin marks, and binding specificities for many TFs determined in vitro . However, such resources are available only for a few major model systems and cannot be easily replicated for new organisms or cell types. Results We present NetProphet 2.0, a ‘data light’ algorithm for TF network mapping, and show that it is more accurate at identifying direct targets of TFs than other, similarly data light algorithms. In particular, it improves on the accuracy of NetProphet 1.0, which used only gene expression data, by exploiting three principles. First, combining multiple approaches to network mapping from expression data can improve accuracy relative to the constituent approaches. Second, TFs with similar DNA binding domains bind similar sets of target genes. Third, even a noisy, preliminary network map can be used to infer DNA binding specificities from promoter sequences and these inferred specificities can be used to further improve the accuracy of the network map. Availability and implementation Source code and comprehensive documentation are freely available at https://github.com/yiming-kang/NetProphet_2.0 . Contact brent@wustl.edu Supplementary information Supplementary dataSupplementary data are available at Bioinformatics online.

Description: Motivation To increase detection power, researchers use gene level analysis methods to aggregate weak marker signals. Due to gene expression controlling biological processes, researchers proposed aggregating signals for expression Quantitative Trait Loci (eQTL). Most gene-level eQTL methods make statistical inferences based on (i) summary statistics from genome-wide association studies (GWAS) and (ii) linkage disequilibrium patterns from a relevant reference panel. While most such tools assume homogeneous cohorts, our G ene-level J oint A nalysis of functional SNPs in C osmopolitan C ohorts (JEPEGMIX) method accommodates cosmopolitan cohorts by using heterogeneous panels. However, JEPGMIX relies on brain eQTLs from older gene expression studies and does not adjust for background enrichment in GWAS signals. Results We propose JEPEGMIX2, an extension of JEPEGMIX. When compared to JPEGMIX, it uses (i) cis-eQTL SNPs from the latest expression studies and (ii) brains specific (sub)tissues and tissues other than brain. JEPEGMIX2 also (i) avoids accumulating averagely enriched polygenic information by adjusting for background enrichment and (ii) to avoid an increase in false positive rates for studies with numerous highly enriched (above the background) genes, it outputs gene q -values based on Holm adjustment of P -values. Availability and implementation https://github.com/Chatzinakos/JEPEGMIX2 . Contact chris.chatzinakos@vcuhealth.org Supplementary information Supplementary dataSupplementary data are available at Bioinformatics online.

Description: Motivation Identification of disease-associated miRNAs (disease miRNAs) is critical for understanding disease etiology and pathogenesis. Since miRNAs exert their functions by regulating the expression of their target mRNAs, several methods based on the target genes were proposed to predict disease miRNA candidates. They achieved only limited success as they all suffered from the high false-positive rate of target prediction results. Alternatively, other prediction methods were based on the observation that miRNAs with similar functions tend to be associated with similar diseases and vice versa. The methods exploited the information about miRNAs and diseases, including the functional similarities between miRNAs, the similarities between diseases, and the associations between miRNAs and diseases. However, how to integrate the multiple kinds of information completely and consider the biological characteristic of disease miRNAs is a challenging problem. Results We constructed a bilayer network to represent the complex relationships among miRNAs, among diseases and between miRNAs and diseases. We proposed a non-negative matrix factorization based method to rank, so as to predict, the disease miRNA candidates. The method integrated the miRNA functional similarity, the disease similarity and the miRNA-disease associations seamlessly, which exploited the complex relationships within the bilayer network and the consensus relationship between multiple kinds of information. Considering the correlation between the candidates related to various diseases, it predicted their respective candidates for all the diseases simultaneously. In addition, the sparseness characteristic of disease miRNAs was introduced to generate more reliable prediction model that excludes those noisy candidates. The results on 15 common diseases showed a superior performance of the new method for not only well-characterized diseases but also new ones. A detailed case study on breast neoplasms, colorectal neoplasms, lung neoplasms and 32 other diseases demonstrated the ability of the method for discovering potential disease miRNAs. Availability and implementation The web service for the new method and the list of predicted candidates for all the diseases are available at http://www.bioinfolab.top . Contact xuanping@hlju.edu.cn or zhang@hlju.edu.cn or lijzh@hit.edu.cn Supplementary information Supplementary dataSupplementary data are available at Bioinformatics online.

Description: Summary The Quest for Orthologs (QfO) is an open collaboration framework for experts in comparative phylogenomics and related research areas who have an interest in highly accurate orthology predictions and their applications. We here report highlights and discussion points from the QfO meeting 2015 held in Barcelona. Achievements in recent years have established a basis to support developments for improved orthology prediction and to explore new approaches. Central to the QfO effort is proper benchmarking of methods and services, as well as design of standardized datasets and standardized formats to allow sharing and comparison of results. Simultaneously, analysis pipelines have been improved, evaluated and adapted to handle large datasets. All this would not have occurred without the long-term collaboration of Consortium members. Meeting regularly to review and coordinate complementary activities from a broad spectrum of innovative researchers clearly benefits the community. Highlights of the meeting include addressing sources of and legitimacy of disagreements between orthology calls, the context dependency of orthology definitions, special challenges encountered when analyzing very anciently rooted orthologies, orthology in the light of whole-genome duplications, and the concept of orthologous versus paralogous relationships at different levels, including domain-level orthology. Furthermore, particular needs for different applications (e.g. plant genomics, ancient gene families and others) and the infrastructure for making orthology inferences available (e.g. interfaces with model organism databases) were discussed, with several ongoing efforts that are expected to be reported on during the upcoming 2017 QfO meeting. Contact selewis@lbl.gov or c.dessimoz@ucl.ac.uk

Description: Motivation In recent years, the massively parallel cDNA sequencing (RNA-Seq) technologies have become a powerful tool to provide high resolution measurement of expression and high sensitivity in detecting low abundance transcripts. However, RNA-seq data requires a huge amount of computational efforts. The very fundamental and critical step is to align each sequence fragment against the reference genome. Various de novo spliced RNA aligners have been developed in recent years. Though these aligners can handle spliced alignment and detect splice junctions, some challenges still remain to be solved. With the advances in sequencing technologies and the ongoing collection of sequencing data in the ENCODE project, more efficient alignment algorithms are highly demanded. Most read mappers follow the conventional seed-and-extend strategy to deal with inexact matches for sequence alignment. However, the extension is much more time consuming than the seeding step. Results We proposed a novel RNA-seq de novo mapping algorithm, call DART, which adopts a partitioning strategy to avoid the extension step. The experiment results on synthetic datasets and real NGS datasets showed that DART is a highly efficient aligner that yields the highest or comparable sensitivity and accuracy compared to most state-of-the-art aligners, and more importantly, it spends the least amount of time among the selected aligners. Availability and implementation https://github.com/hsinnan75/DART Contact hsu@iis.sinica.edu.tw Supplementary information Supplementary dataSupplementary data are available at Bioinformatics online.

Description: Summary Gene expression alterations and potentially underlying gene copy number mutations can be measured routinely in the wet lab, but it is still extremely challenging to quantify impacts of altered genes on clinically relevant characteristics to predict putative driver genes. We developed the R package regNet that utilizes gene expression and copy number data to learn regulatory networks for the quantification of potential impacts of individual gene expression alterations on user-defined target genes via network propagation. We demonstrate the value of regNet by identifying putative major regulators that distinguish pilocytic from diffuse astrocytomas and by predicting putative impacts of glioblastoma-specific gene copy number alterations on cell cycle pathway genes and patient survival. Availability and implementation regNet is available for download at https://github.com/seifemi/regNet under GNU GPL-3 . Contact michael.seifert@tu-dresden.de Supplementary information Supplementary dataSupplementary data are available at Bioinformatics online.

Description: Motivation MicroRNAs (miRNAs) play crucial roles in post-transcriptional regulations and various cellular processes. The identification of disease-related miRNAs provides great insights into the underlying pathogenesis of diseases at a system level. However, most existing computational approaches are biased towards known miRNA-disease associations, which is inappropriate for those new diseases or miRNAs without any known association information. Results In this study, we propose a new method with graph regularized non-negative matrix factorization in heterogeneous omics data, called GRNMF, to discover potential associations between miRNAs and diseases, especially for new diseases and miRNAs or those diseases and miRNAs with sparse known associations. First, we integrate the disease semantic information and miRNA functional information to estimate disease similarity and miRNA similarity, respectively. Considering that there is no available interaction observed for new diseases or miRNAs, a preprocessing step is developed to construct the interaction score profiles that will assist in prediction. Next, a graph regularized non-negative matrix factorization framework is utilized to simultaneously identify potential associations for all diseases. The results indicated that our proposed method can effectively prioritize disease-associated miRNAs with higher accuracy compared with other recent approaches. Moreover, case studies also demonstrated the effectiveness of GRNMF to infer unknown miRNA-disease associations for those novel diseases and miRNAs. Availability and implementation The code of GRNMF is freely available at https://github.com/XIAO-HN/GRNMF/ . Contact luojiawei@hnu.edu.cn Supplementary information Supplementary dataSupplementary data are available at Bioinformatics online.

Description: 〈span〉〈div〉Summary〈/div〉Full-waveform inversion (FWI) of surface waves can provide a valuable contribution to near-surface investigations, since they are mainly sensitive to the S-wave velocity and hold a high signal-to-noise ratio. We investigate the performance of the individual 2-D elastic FWI of Rayleigh and Love waves as well as the feasibility of a simultaneous joint FWI of both wave types. We apply these three methods to synthetic data as well as to field data. In synthetic reconstruction tests we compare the performance of the individual wave type inversions and explore the benefits of a simultaneous joint inversion. In these tests both individual wave type inversions perform similarly well, given that the initial P-wave velocity model is accurate enough. In the case of an inaccurate initial P-wave velocity model, we observe artifacts in the reconstructed S-wave velocity models of the Rayleigh wave FWI and the joint FWI. If the P-wave velocity is sufficiently known the joint FWI can further improve the model of the S-wave velocity. In the field data application the Love wave FWI is superior to the Rayleigh wave FWI, possibly due to the insufficient initial P-wave velocity model, while the joint FWI further improves the inversion result also in this case.〈/span〉

Description: 〈span〉〈div〉Summary〈/div〉As soon as an earthquake starts, the rupture and the propagation of seismic waves redistribute masses within the Earth. This mass redistribution generates in turn a long-range perturbation of the Earth gravitational field, which can be recorded before the arrival of the direct seismic waves. The recent first observations of such early signals motivate the use of the normal mode theory to model the elastogravity perturbations recorded by a ground-coupled seismometer or gravimeter. Complete modeling by normal mode summation is challenging due to the very large difference in amplitude between the prompt elastogravity signals and the direct P-wave signal. We overcome this problem by introducing a two-step simulation approach. The normal mode approach enables a fast computation of elastogravity signals in layered self-gravitating Earth models. The fast and accurate computation of gravity perturbations indicates instrument locations where signal detection may be achieved, and may prove useful in the implementation of a gravity-based earthquake early warning system.〈/span〉

Description: 〈span〉〈div〉Summary〈/div〉Isostatic gravity anomalies provide a measure of the Earth's gravity field free from the gravitational attractions of the topography and its isostatic compensation, most commonly represented by a variation in the depth of a compensating density contrast, for example the Moho. They are used by both geodesists and geophysicists alike, though often for different purposes. Unfortunately though, the effect of subsurface loading on the lithosphere renders transfer function (admittance) methods unusable when surface and subsurface loads coexist. Where they exist, subsurface loads are often expressed in the Bouguer anomaly but not in the topography, and it is shown here that this phase disconnect cannot be faithfully represented by either real- or complex-valued analytic admittance functions. Additionally, many studies that employ the isostatic anomaly ignore the effects of the flexural rigidity of the lithosphere, most often represented as an effective elastic thickness (〈span〉Te〈/span〉), and assume only Airy isostasy, i.e. surface loading of a plate with zero elastic thickness. The consequences of such an omission are studied here, finding that failure to account for flexural rigidity and subsurface loading can result in (1) over- or underestimates of both inverted Moho depths and dynamic topography amplitude, and (2) underestimates of the size of topographic load that can be supported by the plate without flexure. An example of the latter is shown over Europe. Finally, it is demonstrated how low values of the isostatic anomaly variance can actually be biased by these anomalies having low power at the long wavelengths while still possessing high power at middle to short wavelengths, compared to the corresponding Bouguer anomaly power spectrum. This will influence the choice of best-fitting isostatic model if the model is chosen by minimisation of the isostatic anomaly standard deviation.〈/span〉

Description: 〈span〉〈div〉Summary〈/div〉Anisotropic parameters, when considered in a microseismic processing, are typically inverted using perforation shot data and/or simultaneously inverted with microseismic event locations. Inverting using microseismic data alone usually leads to an under-constrained inverse problem that is highly dependent on prior/initial information. We carefully processed the waveforms from perforation shots, and picked P-, SH-, and SV-wave arrival times to mitigate this issue. Because the perforation shot locations are known, the inversion is better-constrained by reducing the number of model parameters while increasing the number of observations. We applied both Maximum A Posteriori (MAP) estimation and Randomized Maximum Likelihood (RML) methods for anisotropic parameter estimation, uncertainty quantification, and trade-off analysis. Results verified the stability of the inversion and revealed the uncertainty and trade-offs among model parameter. In addition, attenuation is generally not considered in microseismic modeling and processing. Our study found that hydraulic stimulation may lead strong increases in seismic attenuation to reservoirs. The attenuation can dramatically change waveform characteristics and cause velocity dispersion. Thus, sonic logs, which are acquired at frequencies much higher than seismic data frequency should not be used directly for data processing in hydraulic stimulated zones.〈/span〉

Description: 〈span〉〈div〉Summary〈/div〉On 26 May 2006 at 23:54 UTC, a moderate shallow crustal earthquake with a moment magnitude of 6.3 occurred in the southern part of Yogyakarta in Java, Indonesia. The earthquake caused severe damages in the area in addition to over 5700 casualties. The cause of this earthquake was initially believed to have been a rupture on the northeast-southwest trending Opak Fault; however, the role of this fault in the earthquake continues to be debated. Therefore, this study presents a subsurface model constructed to characterize the fault geometry associated with the earthquake. We used previously reported aftershock data to image subsurface velocity variations through seismic tomographic inversion of primary waves, shear waves, and their velocity ratio (V〈sub〉p〈/sub〉/V〈sub〉s〈/sub〉). Using data from 10 stations around the hypothetical fault, 588 aftershock events were mostly located 10–15 km east of the Opak River Fault with a maximum depth of approximately 20 km. The seismic tomographic inversion results indicated that severe damage during the earthquake occurred in areas with larger V〈sub〉p〈/sub〉/V〈sub〉s〈/sub〉 ratios associated with unconsolidated sediments, in accordance with previous findings. Furthermore, the configuration of an unnamed fault that was activated during the earthquake is delineated by a velocity anomaly with a depth of up to 5–7 km. This structure is interpreted as a strike slip fault with a reverse component dipping to the east, striking northeast-southwest.〈/span〉

Description: 〈span〉〈div〉Summary〈/div〉Major faults release tectonic stress not only during large earthquakes, but also through slow aseismic slip associated with microseismic activity. Micro earthquakes could be interpreted as the failure of fault asperities. If the slip accumulated during microearthquakes could be estimated with seismological data, the aseismic component of fault slip is more difficult to evaluate. In this study, we are interested in quantifying the aseismic fault slip from the seismic moment released by the simultaneous on fault seismic activity. For that, we study systematically the relationship between total slip and microseismic moment released on a fault modeled as a heterogeneous planar Dieterich-Ruina frictional interface producing earthquakes and aseismic slip. It is shown that for a long enough earthquake sequence, typically involving several ocurrences of the largest asperity failure, the seismic moment released is directly proportional to the total (seismic plus aseismic) slip on the fault accumulated during the sequence. In particular, the total seismic moment divided by the fault area represents a fraction of the total slip that is of the order of the apparent density of asperities. The apparent density of asperities is defined as the fraction of the fault area that experienced a dynamic rupture during the sequence. This result provides a strong mechanical support to extend the use of microseismic activity as creepmeters along major active faults, beyond the common analysis of repeating earthquakes.〈/span〉

Description: 〈span〉〈div〉Summary〈/div〉We introduce in the on-site earthquake early warning (EEW) a partially non-ergodic perspective from the site effects point of view. We consider the on-site EEW approach where the peak ground velocity (PGV) for S-waves is predicted from an early estimate, over the P-waves, of either the peak-displacement (PD) or cumulative squared velocity (IV2). The empirical PD-PGV and IV2-PGV relationships are developed by applying a mixed-effect regression where the site-specific modifications of ground shaking are treated as random effects. We considered a large data set composed of almost 31000 selected recordings in central Italy, a region struck by four earthquakes with magnitude between 6 and 6.5 since the 2009 L’Aquila earthquake. We split the data set into three subsets used for calibrating and validating the on-site EEW models, and for exemplifying their application to stations installed after the calibration phase. We show that the partially non-ergodic models improve the accuracy of the PGV predictions with respect to ergodic models derived for other regions of the world. Moreover, considering PD and accounting for site effects, we reduce the (apparent) aleatory variability of the logarithm of PGV from 0.31-0.36, typical values for ergodic on-site EEW models, to about 0.25. Interestingly, a lower variability of 0.15 is obtained by considering IV2 as proxy, which suggests further consideration of this parameter for the design of on-site EEW systems. Since being site-specific is an inherent characteristic of on-site EEW applications, the improved accuracy and precision of the PGV predicted for a target protection translate in a better customization of the alert protocols for automatic actions.〈/span〉

Description: 〈span〉〈div〉Summary〈/div〉We calculate the change in effusion rate of lava from a volcanic fissure due to pressure changes in the volcanic conduit. The conduit is modelled as a cylinder with elliptical cross section, embedded in an elastic medium. The elliptical shape can represent a wide range of cross sections, according to the value of eccentricity, from almost circular vents to very long and narrow fissures. A 2-D problem is considered assuming invariance of pressure changes and conduit geometry with depth. The problem is solved analytically and expressions for the displacement and the stress fields in the elastic medium are provided. The displacement of the conduit wall is proportional to the ratio between the pressure change and the rigidity of surrounding rocks. The flow rate is a nonlinear function of the pressure change and increases with increasing pressure, due to the elastic deformation of the conduit wall. We consider flow rate oscillations with periods ranging from several minutes to hours, as are often observed during effusive eruptions. Assuming pressure oscillations with these periods, flow rate oscillations resulting from the elastic deformation of the conduit are calculated. The greatest oscillations in flow rate are obtained for very large values of the conduit eccentricity, corresponding to long and narrow volcanic fissures. For example, if a fissure is 100 m long and 2 m large, a pressure oscillation with an amplitude of 1 MPa yields a maximum displacement of the conduit wall equal to about 6 cm and an amplitude of flow rate oscillations of about 20%.〈/span〉

Description: 〈span〉〈div〉Summary〈/div〉We propose a numerical method to compute the inertial modes of a container with near-spherical geometry based on the fully spectral discretisation of the angular and radial directions using spherical harmonics and Gegenbauer polynomial expansion respectively. This allows to solve simultaneously the Poincaré equation and the no penetration condition as an algebraic polynomial eigenvalue problem. The inertial modes of an exact oblate spheroid are recovered to machine precision using an appropriate set of spheroidal coordinates. We show how other boundaries that deviate slightly from a sphere can be accommodated for with the technique of 〈span〉equivalent spherical boundary〈/span〉 and we demonstrate the convergence properties of this approach for the triaxial ellipsoid.〈/span〉

Description: 〈span〉〈div〉Summary〈/div〉We outline a method using Gradient Flow Independent Component Analysis (ICA) to separate signals comprising the coda in a topographically complex setting. We also identify the sources of scattered signals by tracking signal backazimuths over time. The Gradient Flow ICA method is shown to effectively separate signals in the acoustic coda. The method correctly identifies the backazimuth of the first arrival from two 800 kg TNT equivalent explosions as well as subsequent signals scattered by the surrounding topography. Circular statistics is used to analyze the variance, mean, and uniformity of calculated backazimuths. These results have strong implications for understanding the acoustic wavefield by identifying scatterers and inverting for atmospheric conditions.〈/span〉

Description: 〈span〉〈div〉Summary〈/div〉A detailed and precise knowledge of ocean bottom topography is essential in many geoscientific and oceanographic applications. Shipborne echo sounding provides the only direct bathymetric method. However, even after decades of applying this technique only a fraction of the global ocean could be covered. Alternatively, gravity data inversion is a feasible method to infer ocean bottom topography since the gravity field correlates with topography at short to medium wavelengths. Gravity field observables are globally provided by dedicated satellite missions like GOCE, GRACE and GRACE-FO and, over the oceans, by satellite altimetry. Regional and local measurements are realized by means of ground-based, shipborne and airborne gravimetry. For the first time in Europe a jet aircraft was used for airborne gravimetry. During the GEOHALO flight campaign over Italy and the Tyrrhenian, Ionian and Adriatic Seas in June 2012, the German research aircraft HALO carried an entire suite of geodetic-geophysical instrumentation, including gravity meters. The careful processing of the gravity data acquired by the CHEKAN-AM instrument of the German Research Centre of Geosciences allowed to achieve an accuracy at the mGal level. Subsequently, the Parker-Oldenburg inversion was applied to predict ocean bottom topography along the GEOHALO profiles flown over the ocean. To constrain the parameter space in the inversion we used the Crust1.0 model. Finally, the obtained results were compared to the General Bathymetric Chart of the Oceans (GEBCO) in order to estimate the performance of the method. Our study demonstrates that airborne gravimetry aboard a jet aircraft is capable to provide valuable data for regional geoscientific studies.〈/span〉

Description: 〈span〉〈div〉Summary〈/div〉A new method is proposed to estimate separately 3D heterogeneities of intrinsic and scattering attenuation parameters of seismic shear waves. The method, which assumes isotropic scattering, consists of two steps. First, the source, site, and mean path (both intrinsic and scattering attenuation) terms are estimated by an envelope-fitting method for each earthquake event; second, these attenuation parameters are mapped to 3D space in a tomographic manner by using sensitivity kernels of the attenuation parameters. Monte Carlo simulation is adopted to take account of the effect of the depth-dependent velocity structure in the envelope calculation. The proposed method was then applied in south-western Japan. The result showed that both intrinsic and scattering attenuation were stronger in the Kyushu region, where tectonic activity was higher, than in other studied regions. This association of attenuation parameters with tectonic activity is consistent with previous findings. The proposed method thus provides an alternative to multiple lapse time window analysis, a widely used method for separately estimating intrinsic and scattering attenuation factors.〈/span〉

Description: 〈span〉〈div〉Summary〈/div〉We discuss a regularization approach to improving numerical stability of the frequency domain Maxwell equations of electromagnetic geophysics. To enforce divergence-free conditions we add a scaled grad-div operator to the curl-curl equation for electric fields. This deflates the null space of the curl-curl operator and significantly reduces the condition number of the linear system derived for the finite difference (FD) approximation, resulting in faster and more stable iterative solution. We explicitly discuss extensions needed to solve the adjoint problems and consider two apparently different approaches to the sensitivity calculation, which we show are ultimately equivalent. To complete assessment of the new approach in practice we have implemented the modified solver in the ModEM inversion code, and compared it to more standard methods (based on solving with a divergence correction) on inversion of a real data set. This, and simpler tests of the forward solver in isolation, demonstrate that, with appropriate scaling of the added grad-div term, the regularization approach can dramatically improve the efficiency of Krylov subspace methods. Run times for the inversion test are reduced by almost a factor of three, making 3D FD inversion significantly more practical.〈/span〉

Description: 〈span〉〈div〉Summary〈/div〉Anisotropic viscosity is often discussed in the context of flow-induced crystallographic alignment, also known as lattice-preferred orientation, believed to exist throughout the mantle. However, other geophysical systems can produce bulk anisotropic rheology. One example is so-called shape-preferred orientation, such as that imparted by repeated, non-randomly oriented igneous intrusions. Dikes and sills tend to intrude perpendicular to least compressive tectonic stress, and therefore are typically aligned with other coeval intrusions. Intrusions provide distinct planes of weakness, producing, in a bulk sense, an anisotropic fabric to the rocks in which they are intruded. We parameterize this style of anisotropic viscosity with an effective-medium approximation, based on analytically derived characteristics of layered media. We find that the magnitude of anisotropy is highly sensitive to the presence of weak layers, even if they are relatively thin and/or sparsely distributed. We test the effects of our constitutive law on Rayleigh-Taylor instability growth, representing lithospheric instabilities. Our results indicate dramatically decreased stability of lithosphere intruded by dikes and/or sills.〈/span〉

Description: 〈span〉〈div〉Summary〈/div〉We use one month of continuous seismic waveforms from a very dense seismic network to image with unprecedented resolution the shallow damage structure of the San Jacinto fault zone across the Clark fault strand. After calculating noise correlations, high apparent-velocity arrivals coming from below the array are removed using a frequency-wavenumber filter. This is followed by a double-beamforming analysis on multiple pairs of subarrays to extract phase and group velocity information across the study area. The phase and group velocity dispersion curves are regionalized into phase and group velocity maps at different frequencies, and these maps are inverted using a neighbourhood algorithm to build a 3D shear wave velocity model around the Clark fault down to ∼500 m depth. The model reveals strong lateral variations across the fault strike with pronounced low velocity zones corresponding to a local sedimentary basin and a fault zone trapping structure. The results complement previous earthquake- and seismic noise-based imaging of the fault zone at greater depth and clarify properties of structural features near the surface.〈/span〉

Description: 〈span〉〈div〉Summary〈/div〉This paper studies the effect of omission error on height system unification by analyzing the contributions of the low and high spherical harmonic degree components of a local vertical datum (LVD). Previous studies have combined the satellite-only GOCE (Gravity field and steady-state Ocean Circulation Explorer) global geopotential models (GGMs) with a high-resolution GGM (and in some cases topography-implied gravity field signals) in order to account for (or minimize) the omission error of GOCE. This study attempts to answer the question of whether the omission error should be accounted for and to what extent by decomposing the Brazilian Vertical Datum referenced to a tide gauge station at Imbituba (BVD-I) into low and high degree components using the Gaussian averaging function. A combined EIGEN-GOCE GGM is used to compute the high degree datum offset. Overall, the results confirm that (1) the low degree component accounts for the majority of the discrepancy between the BVD-I and the global datum but is insufficient to accurately determine the BVD-I offset parameter, (2) the high degree component contributes little to the overall discrepancy between the datums, but significantly improves the errors in the low degree component, and (3) the signals beyond the resolution of the high-resolution GGM do not significantly impact the datum offset parameters, so that it is not necessary to account for the omission error of the high-resolution GGM using topography-implied gravity field signals. The approach used in this work adds little computational cost, especially over large spatial extents and dense GPS/leveling data, compared to the computation of topography-generated gravity field quantities and has the added advantage of being able to determine the low degree systematic errors in the GPS/leveling data (LVD) using the unbiased satellite-only GOCE models as reference.〈/span〉

Description: 〈span〉〈div〉Summary〈/div〉Earth rotation studies require accurate knowledge of the global oceanic velocity and mass fields, for proper accounting of ocean angular momentum (OAM) effects on the planetary budget. We analyze a new OAM series (1992–2015) based on the solution of a global general circulation model constrained to most existing ocean data. The impact of the data-constrained optimization on OAM is substantial, and particularly essential for calculating effects of global mean ocean mass changes, which can be important for determining annual cycles and long term trends in OAM. The contributions of sea ice to OAM variations, also estimated, are found to be negligible. Uncertainties in OAM series are assessed by comparison with other available estimates. Results indicate low signal-to-noise ratios for all the analyzed OAM series. Comparisons with geodetic, atmospheric and hydrologic data, in the context of the planetary angular momentum budget, point to the continued need for improvements in some or all of the series. Possible paths are offered for producing better OAM estimates in the future.〈/span〉

Description: 〈span〉〈div〉Summary〈/div〉We present 96 new seafloor heat flow determinations, made with a 3.5-m violin-bow probe and collocated with seismic reflection profiles, from the northern Hikurangi margin on the east coast of New Zealand's North Island. Here the Hikurangi Plateau on the Pacific Plate is subducting under the Australian Plate. The background heat flow is 58 ± 8 mW m〈sup〉−2〈/sup〉, consistent with and within the variability of globally observed heat flow (56 ± 15 mW m〈sup〉−2〈/sup〉) for oceanic crust 90 to 120 Ma, the age of the Hikurangi Plateau. Seaward of the deformation front, we find evidence for advective fluid flow associated with basement relief. Landward of the deformation front, we use a two-dimensional steady-state finite element model to quantify the thermal regime. Despite corrections for the effects of bottom water temperature change, bathymetry, and sedimentation there is considerable scatter in the heat flow data including a local and sharp increase in heat flow of up to 35 mW m〈sup〉−2〈/sup〉 observed over the outermost wedge. Variability in the heat flow data is likely due to complex, unmodeled three-dimensional fluid flow. We augment our heat flow measurements with estimates from a bottom-simulating reflection (BSR) and continental bottom hole temperatures and conclude that the effective coefficient of friction,  μ*, is approximately 0.06 in the region of observed slow slip events and increases to 0.18 approximately 50 km landward of the deformation front. This transition in μ* may be marking the downdip edge of overpressures along the subduction thrust, suggesting that slow slip is enabled by overpressure.〈/span〉

Description: 〈span〉〈div〉Summary〈/div〉We analyze noise characteristics of spatial and temporal correlation of 260 continuous GPS sites from Crustal Movement Observation Network of China (CMONOC). These datasets were mainly collected between 2010 and 2016, with an average of 6 years of position time series. In the functional analysis, a clear regional dependence of seasonal movements has been observed and other significant periodic signals are detected nearby the GPS draconitic period and its harmonics. The distribution of these periodic signals shows a spatial correlation, along with non-negligible local inconsistencies. In the stochastic analysis, impacts of the periodicities on the noise assessment have been investigated and Maximum likelihood estimation is used to study noise properties for deseasonalized residual time series having the seasonal signals removed and filtered residual time series having other periodic signals removed further. We demonstrate that for both solutions, the flicker noise is thought to be the dominant time-correlated noise and velocity uncertainties may be underestimated 8–10 times if assuming a pure white noise. Ignoring the periodicities could bias the estimation of noise amplitude, spectral index and velocity uncertainty. After removing the periodic signals, the median flicker noise magnitude shows an average of ∼10 per cent reductions and the noise process shifts closer toward white noise. A correlation between the index variations and RMS variations has been observed, indicating that the index varies more significantly for sites with more periodic signals removed. Besides, the spectral index in the vertical component has a better spatial correlation than that in the horizontal and the spatial distribution of the index of deseasonalized solutions seems to correlate well with the amplitudes of seasonal signals, probably implying common sources of spatial variations of these characteristics. Furthermore, analyses of inter-site correlations indicate that the correlation induced from the periodic signals displays a similar pattern to the deseasonalized solutions, confirming that the period pattern is spatially correlated and can induce time-series correlations. Finally, the stochastic processes of the common mode noise are predominantly featured by spatially correlated flicker and white noise over a wide range, consistent with previous results.〈/span〉

Description: 〈span〉〈div〉Summary〈/div〉Waveform backprojection is a key technique of earthquake-source imaging, which has been widely used for extracting information of earthquake source evolution that cannot be obtained by kinematic source inversion. The technique enjoys considerable popularity, owing to the simplicity of its implementation and the robustness of its processing, but the physical meaning of backprojection images has remained elusive. In this study, we reviewed the mathematical representation of backprojection (BP) and hybrid backprojection (HBP) methods, following the pioneering work of Fukahata et al. (Geophys. J. Int. (2014) 196, 552559), to clarify the physical implications of BP images. We found that signal intensity in BP and HBP images is scaled with the amplitude of the Green’s function that corresponds to a unit-step slip, which results in the signal intensity being depth dependent. We propose variants of BP and HBP, which we call kinematic BP and HBP, respectively, to relate the BP signal intensity to slip motion of an earthquake by modifying the normalizing factors used in the original BP and HBP methods. The original BP and HBP images remain useful for assessing the spatiotemporal strength of the wave radiation, which scales with the amplitude of the Green’s function, whereas the kinematic BP and HBP methods are suitable for imaging the slip motion that is responsible for the high-frequency radiation produced during the source-rupture process.〈/span〉

Description: 〈span〉〈div〉Summary〈/div〉In the field of controlled source seismology, the acoustic 3-D Full Waveform Inversion (FWI) technique has become a common tool for imaging geologically complex structures in land, as well as in marine settings. However, the Earth behaves elastically and therefore, excluding the elastic effect could have a significant impact on the inversion results, especially for interfaces with large S-wave velocity contrasts. To examine the contribution of the elastic approach, we compare acoustic and elastic 3-D FWI applied to a 3-D seismic dataset from the East Pacific Rise (EPR) 9°50′ N, where the subsurface is represented by igneous basaltic rocks. To establish an efficient inversion strategy, we first conducted a number of tests, which suggest a simultaneous, multi-parameter inversion as the most efficient approach when inverting signals with frequencies below 7 Hz. The reduction in the total misfit for the elastic case is 10–15 per cent lower than that for the acoustic one, suggesting that the elastic approach explains the observed data better than the acoustic approach. Furthermore, the compressional velocity images of the upper-oceanic crust obtained using the two approaches differ significantly, not only in velocity magnitude but also structurally. We argue that the results obtained from the acoustic modeling are geologically less plausible and suggest the elastic model as a more reliable representation of the upper-oceanic crust.〈/span〉

Description: 〈span〉〈div〉Summary〈/div〉We investigated an Mw∼6.2 earthquake doublet on the border of the USA and Canada using ALOS2 Light-of-Sight displacements and GPS measurements. We selected three L-band ALOS-2 interfergorams with temporal baselines of one year to extract coseismic deformation maps, in which master and slave images were both acquired in July. A subpixel-based alignment and range spectral splitting techniques under the GAMMA InSAR software framework were applied to improve the interferometric coherence and reduce the effects of phase anomalies in two of the three interferometric pairs due to either ionospheric delay or a potential focusing issues in the generation of the ALOS2 SLC data. The updated interferograms convincingly reveal deformation fringe patterns produced by the two earthquakes. We conducted a nonlinear geophysical inversion to estimate the geometric parameters of the earthquakes with the InSAR and GPS measurements. The best fitting model shows that a thrust faulting on a reverse fault and left-lateral strike-slip faulting on a nearly vertical fault with the centroid depths of 6.7 and 8.5 km, respectively, are most likely responsible for the earthquake doublet. The eastern Denali fault (EDF) and Duke River fault are major active faults in the region and the earthquake doublet could be due to reactivation of the part of the two faults system.〈/span〉

Description: 〈span〉〈div〉Summary〈/div〉We present the theory for and applications of Hamiltonian Monte Carlo (HMC) solutions of linear and nonlinear tomographic problems. HMC rests on the construction of an artificial Hamiltonian system where a model is treated as a high-dimensional particle moving along a trajectory in an extended model space. Using derivatives of the forward equations, HMC is able to make long-distance moves from the current towards a new independent model, thereby promoting model independence, while maintaining high acceptance rates. Following a brief introduction to HMC using common geophysical terminology, we study linear (tomographic) problems. Though these may not be the main target of Monte Carlo methods, they provide valuable insight into the geometry and the tuning of HMC, including the design of suitable mass matrices, and the length of Hamiltonian trajectories. This is complemented by a self-contained proof of the HMC algorithm in the Appendix. A series of tomographic/imaging examples is intended to illustrate (i) different variants of HMC, such as constrained and tempered sampling, (ii) the independence of samples produced by the HMC algorithm, and (iii) the effects of tuning on the number of samples required to achieve practically useful convergence. Most importantly, we demonstrate the combination of HMC with adjoint techniques. This allows us to solve a fully nonlinear, probabilistic traveltime tomography with several thousand unknowns on a standard laptop computer, without any need for supercomputing resources.〈/span〉

Description: 〈span〉〈div〉Summary〈/div〉 〈strong〉Moment Tensors of Hydraulically-Induced AEs:〈/strong〉 Hydraulic fracturing is an important technique in the development of enhanced geothermal systems and unconventional resources. Although the fracture modes induced by hydraulic fracturing influence the recovery efficiency of the resources, the current understanding of this relationship is insufficient. In the present study, we considered the acoustic emissions (AEs) induced during hydraulic fracturing under uniaxial loading conditions in the laboratory, and applied a moment tensor analysis by carefully correcting the coupling condition and directivity of AE transducers. Experiments were conducted for two types of Kurokami-jima granite samples: those with a rift plane perpendicular (Type H) or parallel (Type V) to the expected direction of fracture propagation (i.e. along the loading axis). In the experiments, both sample types experienced a significant number of shear, tensile, and compressive events. The dominant fracture mode for Type H samples is found to be tensile events in which the fracture plane is parallel to the loading axis, whereas for Type V samples, shear events are dominant. This difference suggests that the dominant fracture modes induced by hydraulic fracturing are highly dependent on the relationship between the direction of fracture propagation and orientation of pre-existing weak planes.〈/span〉

Description: 〈span〉〈div〉Summary〈/div〉The size of an earthquake can be defined either from the seismic moment (M〈sub〉0〈/sub〉) or in terms of radiated seismic energy (E〈sub〉r〈/sub〉). These two parameters look at the source complexity from different perspectives: M〈sub〉0〈/sub〉 is a static measure of the earthquake size, whereas E〈sub〉r〈/sub〉 is related to the rupture kinematics and dynamics. For practical applications and for dissemination purposes, the logarithms of M〈sub〉0〈/sub〉 and E〈sub〉r〈/sub〉 are used to define the moment magnitude M〈sub〉w〈/sub〉 and the energy magnitude M〈sub〉E〈/sub〉, respectively. The introduction of M〈sub〉w〈/sub〉 and M〈sub〉E〈/sub〉 partially obscure the complementarity of M〈sub〉0〈/sub〉 and E〈sub〉r〈/sub〉. The reason is due to the assumptions needed to define any magnitude scale. For example, in defining M〈sub〉w〈/sub〉, the apparent stress (i.e. the ratio between M〈sub〉0〈/sub〉 and E〈sub〉r〈/sub〉 multiplied by the rigidity) was assumed to be constant, and under this condition, M〈sub〉w〈/sub〉 and M〈sub〉E〈/sub〉 values would only differ by an off-set which, in turn, depends on the average apparent stress of the analysed dataset. In any case, when the apparent stress is variable and, for example, scales with M〈sub〉0〈/sub〉, the value of M〈sub〉E〈/sub〉 derived from M〈sub〉w〈/sub〉 cannot be used to infer E〈sub〉r〈/sub〉.In this study, we investigate the similarities and differences between M〈sub〉w〈/sub〉 and M〈sub〉E〈/sub〉 in connection with the scaling of the source parameters using a dataset of around 4700 earthquakes recorded at both global and regional scales and belonging to four datasets. These cover different geographical areas and extensions and are composed by either natural or induced earthquakes in the magnitude range 1.5 ≤ M〈sub〉w〈/sub〉 ≤ 9.0. Our results show that M〈sub〉E〈/sub〉 is better than M〈sub〉w〈/sub〉 in capturing the high-frequency ground shaking variability whenever the stress drop differs from the reference value adopted to define M〈sub〉w〈/sub〉. We show that M〈sub〉E〈/sub〉 accounts for variations in the rupture processes, introducing systematic event-dependent deviations from the mean regional peak ground motion velocity scaling. Therefore, M〈sub〉E〈/sub〉 might be a valid alternative to M〈sub〉w〈/sub〉 for deriving ground motion prediction equations for seismic hazard studies in areas where strong systematic stress drop scaling with M〈sub〉0〈/sub〉 are found, such as observed for induced earthquakes in geothermal regions. Furthermore, we analyse the different datasets in terms of their cumulative frequency-magnitude (CFM) distribution, considering both M〈sub〉E〈/sub〉 and M〈sub〉w〈/sub〉. We show that the b values from M〈sub〉w〈/sub〉 (〈span〉b〈/span〉〈sub〉Mw〈/sub〉) and M〈sub〉E〈/sub〉 (〈span〉b〈/span〉〈sub〉ME〈/sub〉) can be significantly different when the stress drop shows a systematic scaling relationship with M〈sub〉0〈/sub〉. We found that 〈span〉b〈/span〉〈sub〉ME〈/sub〉 is nearly constant for all datasets, while 〈span〉b〈/span〉〈sub〉Mw〈/sub〉 shows an inverse linear scaling with apparent stress.〈/span〉

Description: 〈span〉〈div〉Summary〈/div〉Changes in the Earth's magnetic field can deeply modify the polar caps and auroral zones, which are the regions of most frequent precipitation of energetic particles. The present field is characterized by a dominant dipole plus weaker multipolar components. The field varies greatly in time, with the most drastic changes being polarity reversals that take place on average every ∼200,000 years. During a polarity transition the field magnitude may diminish to about 10 per cent of its value prior to the reversal due to a decreasing dipolar component and by becoming mostly multipolar in nature. Polar caps depend on the geomagnetic field configuration so changes in their morphology are expected as a consequence of the variation and reversal of this field. We model polar caps’ location by considering a superposition of the internal geomagnetic field and a uniform external field and then following the open field lines to Earth's surface. Polar caps’ location and shape for different magnetic field reversal scenarios are analyzed in the present work. Two polar caps near the present dipole axis intersection with Earth's surface prevail for a dipole decrease to a certain extent, below which the southern hemisphere polar cap moves to mid-latitudes. An axial dipole collapse gives a pair of polar caps both at mid-latitudes of the southern hemisphere, while in a dipole rotation scenario the polar caps reside at the equator. If reversals occur due to an energy cascade from the dipole to higher degrees, more than two polar caps may appear. In our energy cascade scenario four polar caps at various latitudes of both hemispheres prevail. These results indicate that during reversals auroral zones may reach mid and low-latitude regions, and the atmosphere may become more vulnerable to the direct effect of energetic particle precipitation. This vulnerability is particularly striking at the southern hemisphere where reversed flux patches appear on the core-mantle boundary and weak intensity characterizes the present field at Earth's surface.〈/span〉

Description: 〈span〉〈div〉Summary〈/div〉It is strongly debated whether the interface between the lithosphere and underlying asthenosphere is a temperature-dependent rheological transition, as expected in a thermal convection system, or additionally affected by the presence of melts and/or fluids. Previous surface wave studies of Pacific oceanic lithosphere have found that shear velocity and azimuthal anisotropy vary with seafloor crustal age as expected for a thermal control; however radial anisotropy does not. Various thermo-mechanical models have been proposed to explain this disparate behaviour. Nonetheless, it is unclear how robust the surface wave constraints are, and this is what we test in this study. We apply a Bayesian model space search approach to three published Pacific surface-wave dispersion datasets, two phase-velocity and one combined phase- and group-velocity set, and determine various proxies for the depth of the lithosphere-asthenosphere boundary (LAB) and their uncertainties based on the velocity and radial anisotropy model distributions obtained. In their overall character and pattern with age, the velocity models from different datasets are consistent with each other, although they differ in their values of LAB depths. Uncertainties are substantial (as much as 20 km on LAB depths) and the addition of group-velocity data does not reduce them. Radial anisotropy structures differ even in pattern and display no obvious age dependence. However, given the uncertainties, we cannot exclude that radial anisotropy, azimuthal anisotropy, and velocity models actually reflect compatible, age-dependent, LAB depth estimates. The velocity LAB trends are most like those expected for half-space cooling, because velocity differences persist at old ages, below the depth of common plate cooling models. Any direct signature of sub-ridge melt would be too small-scale to be resolved by these data. However, the velocity-increasing effects of dehydration and depletion due to melting below the ridge could explain why LAB proxy depths tend to a minimum of ∼60 km below young ocean floor.〈/span〉

Description: 〈span〉〈div〉Summary〈/div〉Seismic radiation from indigenous sources can be represented by the excess of model stress over actual stress, a second-order tensor field that Backus & Mulcahy (1976) named the 〈span〉stress glut〈/span〉. We prove a new representation theorem that exactly and uniquely decomposes any stress-glut (or strain-glut) density into a set of orthogonal tensor fields of increasing degree $\alpha$, up to six in number, ordered by their first nonzero polynomial moments. The zeroth-degree field ($\alpha \ = \ 0$) is the projection of the stress-glut density onto its zeroth polynomial moment, which defines the seismic moment tensor, Aki seismic moment ${M_0}$, and centroid-moment tensor (CMT) point source. The higher-degree fields describe mechanism complexity—source variability that arises from the spacetime variations in the orientation of the stress glut, which is the main subject of this theoretical study. The representation theorem generalizes the point-source approximation to a sum of multipoles that features the CMT monopole as its leading term. The first-degree field contributes a dipole tensor with an mechanism orthogonal to the CMT, the second-degree field contributes a quadrupole tensor, and so on, up to six orthogonal fields in all. We define the total scalar moment ${M_{\rm{T}}}$ to be the integral of the scalar moment density, and we use the representation theorem to partition this total moment into a sum of fractional moments for each degree $\alpha $. If the faulting is simple enough, ${M_{\rm{T}}} \approx {M_0}$, and the higher-degree terms will be small. When the faulting is more complex, however, ${M_{\rm{T}}} 〉 {M_0}$; the higher-degree fields will contribute more to the radiation, and this contribution will increase with frequency. Application to simple planar faulting shows that out-of-plane variations in slip-vector orientation reduce ${M_0}/{M_{\rm{T}}}$ more than in-plane variations of similar magnitude, consistent with Kagan's (1988) early study of mechanism complexity. We decompose stress-glut realizations from the Graves & Pitarka (2016) rupture simulator; typical values of ${M_0}/{M_{\rm{T}}}$ are 0.82–0.92, consistent with analytical results. The higher-degree fields of the Graves-Pitarka sources typically radiate up to $\alpha \ = \ 4$; only the isotropic term is zero. We compute synthetic seismograms to illustrate the radiation patterns of the higher-degree fields and their frequency dependence. Decomposition of the Duputel & Rivera (2017) source model for the 2016 Kaikoura earthquake (M〈sub〉w〈/sub〉 7.8) indicates that the radiation from the higher-degree fields was large enough (${M_0}/{M_{\rm{T}}}$ = 0.82) that it may be possible to invert global datasets for low-degree multipoles. Other applications of the stress-glut representation theorem are discussed.〈/span〉

Description: 〈span〉〈div〉SUMMARY〈/div〉We present a joint analysis of newly acquired gravity and teleseismic data in the North Tanzanian Divergence, where the lithospheric break-up is at its earliest stage. The impact of a mantle upwelling in more mature branches of the East African Rift has been extensively studied at a lithospheric scale. However, few studies have been completed that relate the deep-seated mantle anomaly detected in broad regional seismic tomography with the surface deformation observed in the thick Archaean Pan-African suture zone located in North Tanzania. Our joint inversion closes the gap between local and regional geophysical studies, providing velocity and density structures from the surface down to 〈span〉ca〈/span〉. 250  km depth with new details. Our results support the idea of a broad mantle upwelling rising up to the lithosphere and creating a thermal modification along its path. However, our study clearly presents an increasing amplitude of the associated anomaly both in velocity and density above 200 km depth, which cannot be solely explained by a temperature rise. We infer from our images the combined impact of melt (2–3 per cent), composition and hydration that accompany the modification of a thick heterogenous cratonic lithosphere are a response to the hot mantle rising. The detailed images we obtained in density and velocity assert that Archaean and Proterozoic units interact with the mantle upwelling to restrict the lithosphere modifications within the Magadi–Natron–Manyara rift arm. The composition and hydration variations associated with those units equilibrate the thermal erosion of the craton root and allow for its stability between 100 and 200  km depth. Above 80 km depth, the crustal part is strongly affected by intruding bodies (melt and gas) which produces large negative anomalies in both velocity and density beneath the main magmatic centres. In addition to the global impact of a superplume, the velocity and density anomaly pattern suggests a 3-D distribution of the crust and mantle lithospheric stretching, which is likely to be controlled by inherited fabrics and enhanced by lateral compositional and hydration variations at the Tanzanian craton-orogenic belt boundary.〈/span〉

Description: 〈span〉〈div〉SUMMARY〈/div〉The nonlinear interaction of ocean surface waves produces coherent infrasound noise—microbaroms—between 0.1 and 0.5 Hz. Microbaroms propagate through the atmosphere over thousands of kilometres due to low absorption and efficient ducting between the ground and the stratopause. These signals are globally and permanently detected by the International Monitoring System (IMS) infrasound network, which has been established to monitor compliance with the Comprehensive Nuclear-Test-Ban Treaty. At the International Data Centre (IDC) in Vienna, where IMS data are routinely processed, microbarom detections appear in overlapping frequency bands, and are treated as false alarms. Therefore, understanding the variability in microbarom detections is essential to support the IDC in the reduction of the false alarm rate. In this study, microbarom amplitudes and the direction of arrivals at the German infrasound station IS26 were modelled. For the simulations, the source was described by an operational ocean wave interaction model, and the signal amplitude was modelled using a semi-empirical attenuation relation. This relation strongly depends on middle atmosphere (MA; i.e. 15–90 km altitude) dynamics; however, vertical temperature and wind profiles, provided by numerical weather prediction (NWP) models, exhibit significant biases and differences when compared with high-resolution light detection and ranging instrument (lidar) soundings in altitudes where infrasound signals propagate. To estimate uncertainties in the modelled amplitude, a fully autonomous lidar for MA temperature measurements was installed at IS26. Temperature and wind perturbations, considering observed biases and deviations, were added to the operational high-resolution atmospheric model analysis produced by the European Centre for Medium-Range Weather Forecasts. Such uncertainties in horizontal winds and temperature strongly impact propagation conditions, explaining almost 97 per cent of the actual detections, compared to 77 per cent when using the direct output of the NWP model only. Incorporating realistic wind and temperature uncertainties in NWP models can thus significantly improve the understanding of microbarom detections as well as the detection capability of a single station throughout the year.〈/span〉

Description: 〈span〉〈div〉SUMMARY〈/div〉In the field of controlled source seismology, the acoustic 3-D Full Waveform Inversion (FWI) technique has become a common tool for imaging geologically complex structures in land, as well as in marine settings. However, the Earth behaves elastically and, therefore, excluding the elastic effect could have a significant impact on the inversion results, especially for interfaces with large 〈span〉S〈/span〉-wave velocity contrasts. To examine the contribution of the elastic approach, we compare acoustic and elastic 3-D FWI applied to a 3-D seismic data set from the East Pacific Rise (EPR) 9°50′ N, where the subsurface is represented by igneous basaltic rocks. To establish an efficient inversion strategy, we first conducted a number of tests, which suggest a simultaneous, multiparameter inversion as the most efficient approach when inverting signals with frequencies below 7 Hz. The reduction in the total misfit for the elastic case is 10-15 per cent lower than that for the acoustic one, suggesting that the elastic approach explains the observed data better than the acoustic approach. Furthermore, the compressional velocity images of the upper-oceanic crust obtained using the two approaches differ significantly, not only in velocity magnitude but also structurally. We argue that the results obtained from the acoustic modelling are geologically less plausible and suggest the elastic model as a more reliable representation of the upper-oceanic crust.〈/span〉

Description: 〈span〉〈div〉SUMMARY〈/div〉We determine 〈span〉m〈/span〉〈sub〉B〈/sub〉, the original body wave magnitude developed by Gutenberg and Richter over the period 1942–1956, for about 3300 〈span〉M〈/span〉〈sub〉w〈/sub〉 ≥ 6 earthquakes for the period 1988–present using modern broad-band seismograms. The main objective is to extend the database of energy-related parameters by combining 〈span〉m〈/span〉〈sub〉B〈/sub〉 databases for recent and old events. The radiated energy 〈span〉E〈/span〉〈sub〉R_B〈/sub〉 (in erg) computed from 〈span〉m〈/span〉〈sub〉B〈/sub〉 using the Gutenberg & Richter relation $\log {E_{\mathrm{ R}\_\mathrm{ B}}} = 2.4{m_\mathrm{ B}} + 5.8$ agrees very well with 〈span〉E〈/span〉〈sub〉R〈/sub〉 estimated with modern techniques, especially for large deep earthquakes. Thus, 〈span〉E〈/span〉〈sub〉R_B〈/sub〉 is useful as a proxy for 〈span〉E〈/span〉〈sub〉R〈/sub〉 to investigate the global diversity of earthquake characteristics and physics over an extended period of time.〈/span〉

Description: 〈span〉〈div〉SUMMARY〈/div〉Locating microseismic events is essential for many areas of seismology including volcano and earthquake monitoring and reservoir engineering. Due to the large number of microseismic events in these settings, an automated seismic location method is required to perform real time seismic monitoring. The measurement environment requires a precise and noise-resistant event location method for seismic monitoring. In this paper, we apply Multichannel Coherency Migration (MCM) to automatically locate microseismic events of induced and volcano-tectonic seismicity using sparse and irregular monitoring arrays. Compared to other migration-based methods, in spite of the often sparse and irregular distribution of the monitoring arrays, the MCM can show better location performance and obtain more consistent location results with the catalogue obtained by manual picking. Our MCM method successfully locates many triggered volcano-tectonic events with local magnitude smaller than 0, which demonstrates its applicability on locating very small earthquakes. Our synthetic event location example at a carbon capture and storage site shows that continuous and coherent drilling noise in industrial settings will pose great challenges for source imaging. However, automatic quality control techniques including filtering in the frequency domain and weighting are used to automatically select high-quality data, and can thus effectively reduce the effects of continuous drilling noise and improve source imaging quality. The location performance of the MCM method for synthetic and real microseismic data sets demonstrates that the MCM method can perform as a reliable and automatic seismic waveform analysis tool to locate microseismic events.〈/span〉

Description: 〈span〉〈div〉Abstract〈/div〉In order to investigate the subduction and continental collision of the North China and Yangtze blocks, two magnetotelluric profiles were obtained across the Dabie Orogenic Belt, the Lower Yangtze Depression, and the Jiangnan Orogenic Belt in the central section (Susong-Anqing section) of the middle-lower Reaches of the Yangtze River in China. After data processing and inversions, we obtained electrical models of the crust and upper mantle. The prominent feature revealed by the inversion is an extensive, arched conductive layer that extends from the middle-lower crust to the upper mantle. To the southeast of this layer, another wedge-shaped conductor is located beneath the Lower Yangtze Depression and the Jiangnan Orogenic Belt. In addition, several resistors, which are distributed from the lower crust to the upper mantle, were also revealed by these lines. These resistors are beneath the conductive layer and separated by vertical conductive bands. Based on these electrical structures, we identified several major faults, including the Tanlu Fault in the eastern part of the Dabie Orogenic Belt, the Yangtze Deep Fault, and the main thrust fault in the Lower Yangtze Depression, which are middle-upper crustal faults. In addition, a ‘Crocodile’ structure was revealed by the major faults in the depression, which are connected by a middle-lower crustal detachment and the surrounding resistive strata. Based on the different electrical structures of the three belts and the results of previous studies, we conclude that subsequent to the slab subducting toward the North China Block during the Triassic, the subduction-collision process that occurred in the study area can be divided into three stages. In the first stage, the weak layer and the Yangtze Fold Belt formed during subduction of the Yangtze Block beneath the North China Block, and the Dabie Orogenic Belt formed during the collision process. In the second stage, the slab buckled as deep material upwelled, and the ‘Crocodile’ structure formed due to the weak layer. In the third stage, the breakup and sinking of the slab caused the Yangtze Fold Belt to subside.〈/span〉

Description: 〈span〉〈div〉Summary〈/div〉It is well known that for a seismic source on or near a material discontinuity, the moment tensor is ambiguous in the sense that its value depends on the exact location of the source. In this short paper we demonstrate that this ambiguity can be resolved by considering components of the source moment and potency tensors separately. For a general source and any location ‘at’ an interface, certain components of its moment tensor and certain components of its potency tensor can be determined unambiguously. The other components depend on the exact location of the source. To investigate these results it is convenient to allow for an asymmetric moment tensor, where the antisymmetric part is physically equivalent to a torque source. These results are illustrated using ray theory for a source at a strong material contrast.〈/span〉

Description: 〈span〉〈div〉SUMMARY〈/div〉Earth's large-scale topography and lower mantle structure are linked to past tectonic motions and mantle flow, making it possible to gain insights in the properties of the solid Earth from time-dependent global convection models driven by tectonic reconstructions. Recent work suggests that the amplitude of residual topography, obtained by subtracting models of isostatic topography from total topography, may be up to ∼1,000 m at spherical harmonic degree two (wavelength ∼13,300 km) and 〉 200 m at spherical harmonic degree 12 (wavelength ∼3,000 km). The amplitude of dynamic topography and the structure of the lower mantle predicted by time-dependent forward mantle flow models both depend on the physical assumption and on model parameters. Here, I investigate the consequence of using the Boussinesq, extended Boussinesq or truncated anelastic liquid approximation (TALA) in time-dependent flow models for predicting present-day mantle structure and dynamic topography; I characterise the sensitivity of the spectral distribution of dynamic topography amplitude to the boundary conditions and model setup for the computation of dynamic topography; and I investigate the sensitivity of model results to parameters including the depth- and temperature-dependence of viscosity, the model initial age, and the density of the basal layer. Extended Boussinesq and TALA models are preferred to Boussinesq models that overpredict the volume of lower mantle slabs and the amplitude of long-wavelength dynamic topography. The correlation between dynamic topography and residual topography models at degrees one-three generally ranges between 0.4–0.5. The flow models better predict the geographical location of Large Low Shear Velocity Provinces (LLSVPs) than that of lower mantle slabs, which cover smaller areas in map view. I show that preserving shallow lateral viscosity variations in the computation of dynamic topography increases the amplitude of dynamic topography for wavelengths 〉 6,000 km. Parameter trade-offs exist to fit both deep and surface constraints. The best-fit model cases considered either a moderately dense basal layer (approximately 2 per cent denser than ambient mantle) or weak temperature- and pressure-dependence of lower mantle viscosity. The fit to present-day constraints does not significantly deteriorate when extending the reconstruction of mantle flow from 200 Ma to 410 Ma, reflecting that seismic tomography models capture the history of mantle convection over the last 200–250 Myr, and suggesting that paleogeographically-constrained mantle flow models should be compared to time-dependent surface geological constraints.〈/span〉

Description: 〈span〉〈div〉SUMMARY〈/div〉Guided waves in a water layer overlaying an elastic half-space are known as normal modes. They are often present in seismic recordings at long offsets in shallow-water environment and generally considered coherent noise. The normal modes, however, carry important information about the near-surface and, as demonstrated by a number of authors, can be used to obtain the shallow velocity model. There is a growing evidence that the latter needs not to be isotropic due to various geological reasons. Motivated by that, we consider the normal-mode propagation in case the elastic half-space exhibits orthorhombic anisotropy. We derive the period equation that describes the normal-mode phase velocity dispersion. To simplify the complicated expression, we present acoustic and ellipsoidal orthorhombic approximations. We also outline the approach towards the group velocity and group azimuth calculation and apply it to the ellipsoidal case to obtain concise and intuitive expressions. Using numerical test, we study the relation between phase and group domains in elastic orthorhombic case. The deviation between velocities and azimuths in these domains is the strongest for low frequencies and it rapidly decreases with increasing frequency. For higher frequencies, the anisotropy effects of the underlaying half-space are barely detectable since the observed signal is composed mainly of the direct acoustic wave, resulting in the two domains being nearly indistinguishable.〈/span〉

Description: 〈span〉〈div〉SUMMARY〈/div〉We propose a theory for rotational and strain ambient noise interferometry, motivated by the recent development of rotational ground motion sensors and distributed acoustic sensing (DAS) technology. In this context, we demonstrate that displacement, strain and rotation interferograms can be generically written in the form of a representation theorem, that is, as a solution to the seismic wave equation that we refer to as the interferometric wavefield. The physical quantity (displacement, strain or rotation) determines the distributed source of the interferometric wavefield, as well as an observational operator that extracts the correct type of noise correlation function. The proposed interferometric equations are free of assumptions on the distribution of noise sources or the equipartitioning of the ambient field, typically required for Green’s function retrieval. In addition to being valid for any kind of heterogeneous source and visco-elastic medium, they allow us to account for measurement details, such as the gauge length in DAS. We illustrate the practical feasibility of our approach with a series of numerical examples, based on regional-scale, spectral-element simulations of the interferometric wavefield. Specifically, we compare displacement and strain interferograms for homogeneous and heterogeneous Earth models, and for homogeneous and heterogeneous noise sources. Ultimately, our developments are intended to enable adjoint-based waveform inversion with emerging measurement technologies that provide spatial gradient information in addition to conventional seismic displacement recordings.〈/span〉

Description: 〈span〉〈div〉Summary〈/div〉With a Hamilton-Jacobi equation in Cartesian coordinates as a starting point, it is common to use a system of ordinary differential equations describing the continuation of first-order derivatives of phase-space perturbations along a reference ray. Such derivatives can be exploited for calculation of geometrical spreading on the reference ray, and for establishing a framework for second-order extrapolation of traveltime to points outside the reference ray. The continuation of first-order derivatives of phase-space perturbations has historically been referred to as dynamic ray tracing. The reason for this is its importance in the process of calculating amplitudes along the reference ray. We extend the standard dynamic ray tracing scheme to include higher-order derivatives of the phase-space perturbations. The main motivation is to extrapolate and interpolate amplitude and phase properties of high-frequency Green’s functions to nearby (paraxial) source and receiver locations. Principal amplitude coefficients, geometrical spreading factors, geometrical spreading matrices, ray propagator matrices, traveltimes, slowness vectors, and curvature matrices are examples of quantities for which we enhance the computation potential. This, in turn, has immediate applications in modelling, mapping, and imaging. Numerical tests for 3D isotropic and anisotropic heterogeneous models yield clearly improved extrapolation results for traveltime and geometrical spreading. One important conclusion is that the extrapolation function for geometrical spreading must be at least third order to be appropriate at large distances away from the reference ray.〈/span〉

Description: 〈span〉〈div〉Summary〈/div〉Of the several factors involved in the development of magnetic fabrics in fault zones at shallow crustal levels, lithology and deformation intensity have probably the most important consequences for the reconstruction of their kinematic history. The basement-involved Cenozoic thrusts in the Demanda Massif (N Spain) provide the opportunity for testing the applicability of Anisotropy of Magnetic Susceptibility (AMS) to the study of deformation in cataclastic fault rocks belonging to shallow fault zones. The Rastraculos thrust is a relatively minor basement thrust (dip-slip movement of 2 km defined from cross-sections and geological maps) of Cenozoic age. This thrust contains a re-activated fault zone involving different rock types both belonging to its hangingwall (Paleozoic) and its footwall (Triassic sandstones and dolostones and Jurassic limestones). AMS results show magnetic foliations parallel or slightly oblique to the fault zone, and both transport-parallel (projected onto the foliation plane) and transport-perpendicular (parallel to the observed intersection lineation) magnetic lineations. The two types of strain/magnetic fabric relationships can be related to deformational and mineralogical features inferred from the direct analysis of thin and polished sections under the microscope and the naked eye, respectively. Analysis of fault rocks in the Rastraculos fault zone indicates that in cataclasites, magnetic fabrics are particularly dependent on lithology and hence magnetic mineralogy. The results obtained prove the usefulness of AMS in fault zones where kinematic indicators are scarce and also give clues on the number of samples necessary to define magnetic susceptibility axes, depending on grain size, ellipsoid shapes and magnetic mineralogy.〈/span〉

Description: 〈span〉〈div〉Summary〈/div〉The subsurface offset is a possible link connecting wave-equation migration methods, such as reverse-time migration, with the angle-domain. However, it describes an artificial reflection geometry that splits the reflection point between incident and scattered waves. The split configuration contradicts basic principles of continuum mechanics, as it represents action-at-a-distance, and has kinematics distinct from those of ordinary physical reflection. The angle-domain image, computed as a Radon transform of the subsurface offset image, implicitly inherits the split configuration, which distinguishes it from conventional angle-domain decomposition. As a consequence, it shows a dissimilar moveout response to erroneous migration velocity, but still carries valuable information about the migration velocity error.Conventional linearized traveltime inversion techniques, applied for velocity optimization, are most likely to fail while inverting angle-domain images produced from subsurface offset split reflections. These images should be labeled differently to highlight the split mechanism, and input to traveltime inversion algorithms only after a proper generalization. We present a modified linearized traveltime inversion formulation, appropriate for split reflections in the image space. The key ingredient in the algorithm is a subsurface offset depended extension of the traveltime change along migrated split reflection rays. A generalized reflection tomography scheme is proposed accordingly to tie imaging errors, in either the subsurface offset domain or its equivalent angle-domain map, to the migration velocity error.〈/span〉

Description: 〈span〉〈div〉Summary〈/div〉Surface wave inversion based on higher modes of dispersion curves is increasingly being applied to estimate near-surface S-wave velocity structures using active or passive sources. Moreover, sensitivity analysis of dispersion curves is important for estimating the confidence of inverted results. In this study, we define a sensitivity function and a sensitivity kernel, give their analytic expressions, and then propose a method to estimate the confidence of inverted results and redesign the initial structure through sensitivity analysis. Three kinds of common layered models (velocity-increasing-with-depth (regular), high-velocity-layer (HVL), and low-velocity-layer (LVL)) are considered. Results show that dispersion curves are more sensitive to S-wave velocity than P-wave velocity and density. The frequency band of the sensitivity function increases as the modes of dispersion curves increase, and the higher modes can penetrate additional depth. For the HVL model, the dispersion curves are particularly sensitive to those of layers above the high-velocity layer. By contrast, for the LVL model, the dispersion curves are particularly sensitive to the low-velocity layer. Aside from the analysis of sensitivity, the sensitivity functions can also be used for inversion. Thus, we perform an S-wave velocity inversion for synthetic records calculated from the structural model of Shanghai, China. Synthetic noise records are converted into the frequency-velocity domain (c-f panel) by the vector wavenumber transform method (VWTM). Then, the dispersion curves of multi-mode are determined from the peaks of the c-f panel. The inverted structure perfectly matches the target structure except the deeper structures that is consistent with the sensitivity analysis of the determined dispersion curves. Finally, we apply the same procedure to the real case in Linyi, Shandong, China. The first three modes of dispersion curves accurately match the data. The synthetic and real examples indicates that the combined use of higher modes offers additional constraint for the velocity structure, which leads to better results of inversion than that of only the fundamental mode. In addition, the sensitivity analysis offers us a tool to estimate the confidence of the inverted structure and redesign the initial structure.〈/span〉

Description: 〈span〉〈div〉Summary〈/div〉In this paper we present a series of mathematical abstractions for seismologically relevant wave equations discretized using finite-element methods, and demonstrate how these abstractions can be implemented efficiently in computer code. Our motivation is to mitigate the combinatorial complexity present when considering geophysical waveform modeling and inversion, where a variety of spatial discretizations, material models, and boundary conditions must be considered simultaneously. We accomplish this goal by first considering three distinct classes of abstract mathematical models: (1) those representing the physics of an underlying wave equation, (2) those describing the discretization of the chosen equation onto a finite-dimensional basis, and (3) those describing any spatial transforms. A full representation of the discrete wave equation can then be constructed using a hierarchical nesting of models from each class. Additionally, each class is functionally orthogonal to the others, and with certain restrictions models within one class can be interchanged independently from changes in another. We then show how this re-casting of the relevant equations can be implemented concisely in computer software using an abstract object-oriented design, and discuss how recent developments in the numerical and computational sciences can be naturally incorporated. This builds to a set of results where we demonstrate how the developments presented can lead to an implementation capable of multi-physics waveform simulations in completely unstructured domains, on both hypercubical and simplical spectral-element meshes, in both two and three dimensions, while remaining concise, efficient, and maintainable.〈/span〉

Description: 〈span〉〈div〉Summary〈/div〉As temperature increases with depth and the creep resistance of rock decreases exponentially, a high viscosity sub-lithospheric layer, just beneath the “elastic” lithosphere is expected to exist. Depending on the temperature profile, a low viscosity asthenosphere may also exist if the temperature deeper down gets high enough. Since the temperature profile is expected to change laterally - especially from below the oceans to cratonic areas underneath continents, rock properties of the lithosphere, high viscosity sub-lithosphere and low viscosity asthenosphere are expected to change laterally. Our aim is to constrain sub-lithospheric properties (depth, thickness and viscosity), lateral lithospheric thickness variations and asthenospheric properties using observed GIA data. A Coupled Laplace-Finite Element Method is used to compute gravitationally self-consistent sea level with time dependent coastline and rotational feedback in addition to changes in deformation, gravity and the state of stress. We start with the VM5a-ICE-6G_C model combination and then modify the lithospheric, sub-lithospheric and asthenospheric properties (including lateral thickness variation) while keeping the mantle viscosities the same as VM5a. Through this study, we confirm that the sub-lithospheric and asthenospheric properties can significantly affect the predicted global RSL, present-day gravity rate-of-change (g-dot) and uplift rate (u-dot) in Laurentia and Fennoscandia. In addition, incorporating the elastic lithosphere with lateral thickness variation, sub-lithosphere and asthenosphere can improve the fit to global RSL, but the predicted peak values of g-dot and u-dot in Laurentia may decrease slightly but not significant enough to affect the fit to the observed data. Our results prefer an elastic lithosphere that has maximum thickness of 140 km under continental cratons but reduces to 60 km underneath the oceans. The results preferred depth of the asthenospheric bottom is around 190-200 km with asthenospheric viscosity around 10〈sup〉20〈/sup〉 Pa s. Finally, we show that the best laterally heterogeneous mantle model we found in previous publication when combined with the lateral thickness varying lithosphere gives the best fit to global RSL and peak g-dot and u-dot in Laurentia simultaneously.〈/span〉

Description: 〈span〉errata, addendaCore flowData assimilationError estimationStochastic modelsKalman Filtercore dynamicsgeomagnetic field modellingsatellite observations〈/span〉

Description: 〈span〉〈div〉Summary〈/div〉Shales commonly exhibit anisotropy in their elastic wave velocity, which directly impacts the accuracy of seismic imaging and their geomechanical response to drilling and completions. Anisotropy is often caused by mineralogical layering, fractures, and rock fabric (i.e. oriented grains and intrinsic anisotropy of clay sediments). However, the relative impact of each of these features on macroscopic shale properties is not well understood. We combined scanning electron microscopy (SEM) and X-ray micro-computed tomography (CT) to image the mineralogical and structural heterogeneity of Mancos Shale and converted the acquired CT and SEM images into heterogeneous 2D elastic models. We used wave propagation numerical simulations to understand the effects that layering and fractures have on elastic wave velocity anisotropy. Consistent 〈span〉Vp〈/span〉/〈span〉Vs〈/span〉 ratios around 1.46 for modeled and measured velocities validates SEM observations of a quartz-dominated shale lithology. CT-derived models containing layering and fractures exhibit 28.6 per cent and 58.8 per cent of the 〈span〉Vp〈/span〉 and 〈span〉Vs〈/span〉 anisotropy observed in the laboratory, whereas SEM derived models exhibit 74.5 per cent and 73.2 per cent of the anisotropy, respectively. The increased anisotropy of SEM-derived elastic models is a result of the ability of the SEM to discern individual mineral grains and microstructural features, whereas the CT models require the use of an effective medium theory to model variations of lithology. Overall, modeled wave propagation perpendicular to bedding more closely captures the experimental velocities than parallel to bedding. Therefore, sub-resolution rock fabric anisotropy likely accounts for the relatively larger velocity mismatch in the parallel direction, and is likely responsible for the decreased anisotropy in coarse rock models. Future modeling would require higher resolution images to structurally constrain these features and/or anisotropic elements to account for fabric anisotropy. Despite some limitations, our study provides a reliable procedure to estimate anisotropy of dynamic mechanical properties of laminated shales using SEM and CT imaging combined with numerical simulation of wave propagation.〈/span〉

Description: 〈span〉〈div〉Summary〈/div〉We implement an analytical model based on flexural deflection of a thin elastic disc to investigate the magnitude of lithospheric decompression caused by deglaciation at upper crustal magmatic reservoirs. Considering a published numerical climate model describing the space-time evolution of deglaciation after the Last Glacial Maximum (LGM) along the Southern Volcanic Zone (SVZ) of the Andes, we demonstrate that changes in pressure at upper crustal levels (〈 10 km depth) at the scale of several hundred years are of the order of 10–100 MPa. Total decompression and decompression rate (300–150 kPa/yr) are 1–2 orders of magnitude larger than values previously estimated by other authors who assume that glacial loads are supported by an elastic half space, i.e. of infinite elastic thickness. The large decompression caused by flexural unbending of an elastic plate of finite thickness as assumed here can easily surpass the tensile strength of rocks (5–20 MPa), creating adequate conditions for failure of the reservoir walls, dike propagation inside and outside the reservoir and the eventual collapse of the reservoir accompanying an explosive eruption. We apply our results to the analysis of post-glacial eruptions of SVZ volcanoes, which erupted large volumes (〉 10 km〈sup〉3〈/sup〉) of mafic ignimbrites hundreds to one thousand years after deglaciation onset. We show that this time lag is necessary to achieve a decompression of several tens of megapascals at depths of several kilometers that are consistent with the location of magmatic reservoirs as estimated by independent petrologic, seismic and/or geodetic studies. Moreover the northward increase of this time lag is in agreement with a smaller size of the Andean ice cap in the north than in the south during the LGM. For wet, volatile-rich magmas typical of subduction zones, the effect of large decompression at upper crustal reservoirs caused by flexural unbending of the lithosphere after deglaciation could play a major role in promoting large explosive eruptions through devolatization of the magma, during past deglaciation events as demonstrated here for the LGM along the SVZ and current accelerated ice retreat caused by climate change over large segments of subduction-related arcs at higher latitudes.〈/span〉

Description: 〈span〉〈div〉Summary〈/div〉The design of an array configuration is an important task in array seismology during experiment planning. Often the array response function (ARF), which depends on the relative position of array stations and frequency content of the incoming signals, is used as the array design criterion. In practice, additional constraints and parameters have to be taken into account, e.g. land ownership, site-specific noise levels or characteristics of the seismic sources under investigation. In this study a flexible array design framework is introduced which implements a customizable scenario modelling and optimization scheme by making use of synthetic seismograms. Using synthetic seismograms to evaluate array performance makes it possible to consider additional constraints. We suggest to use synthetic array beamforming as an array design criteria instead of ARF. The objective function of the optimization scheme is defined according to the monitoring goals, and may consist of a number of sub-functions. The array design framework is exemplified by designing a 7-station small scale array to monitor earthquake swarm activity in North-West Bohemia / Vogtland in central Europe. Two sub-functions are introduced to verify the accuracy of horizontal slowness estimation; one to suppress aliasing effects due to possible secondary lobes of synthetic array beamforming calculated in horizontal slowness space, and the other to reduce the event’s mislocation caused by miscalculation of the horizontal slowness vector. Subsequently, a weighting technique is applied to combine the sub-functions into one single scalar objective function to use in the optimization process.〈/span〉

Description: 〈span〉〈div〉Summary〈/div〉Traveltime approximation plays an important role in seismic data processing, e. g. anisotropic parameter estimation and seismic imaging. By exploiting seismic traveltimes, it is possible to improve the accuracy of anisotropic parameter estimation and the resolution of seismic imaging. Conventionally, the traveltime approximations in anisotropic media are obtained by expanding the anisotropic eikonal equation in terms of the anisotropic parameters and the elliptically anisotropic eikonal equation based on perturbation theory. Such an expansion assumes a small perturbation and weak anisotropy. In a realistic medium, however, the assumption of small perturbation likely breaks down. We present a retrieved zero-order deformation equation that creates a map from the anisotropic eikonal equation to a linearized partial differential equation system based on the homotopy analysis method. By choosing the linear and nonlinear operators in the retrieved zero-order deformation equation, we develop new traveltime approximations that allow us to compute the traveltimes for a medium of arbitrarily strength anisotropy. A comparison of the traveltimes and their errors from the homotopy analysis method and from the perturbation method suggests that the traveltime approximations provide a more reliable result in strongly anisotropic media.〈/span〉

Description: 〈span〉〈div〉Summary〈/div〉Reservoir models are numerical representations of the subsurface petrophysical properties such as porosity, volume of minerals and fluid saturations. These are often derived from elastic models inferred from seismic inversion in a two-step approach: first, seismic reflection data are inverted for the elastic properties of interest (such as density, P-wave and S-wave velocities); these are then used as constraining properties to model the subsurface petrophysical variables. The sequential approach does not ensure a proper propagation of uncertainty throughout the entire geo-modelling workflow as it does not describe a direct link between the observed seismic data and the resulting petrophysical models. Rock physics models link the two domains. We propose to integrate seismic and rock physics modelling into an iterative geostatistical seismic inversion methodology. The proposed method allows the direct inference of the porosity, volume of shale and fluid saturations by simultaneously integrating well-logs, seismic reflection data and rock physics model predictions. Stochastic sequential simulation is used as the perturbation technique of the model parameter space, a calibrated facies-dependent rock physics model links the elastic and the petrophysical domains and a global optimizer based on cross-over genetic algorithms ensures the convergence of the methodology from iteration to iteration. The method is applied to a 3D volume extracted from a real reservoir dataset of a North Sea reservoir and compared to a geostatistical seismic AVA.〈/span〉

Description: 〈span〉〈div〉Summary〈/div〉We determine 〈span〉mB〈/span〉, the original body-wave magnitude developed by Gutenberg and Richter over the period 1942–1956, for about 3,300 〈span〉Mw 〈/span〉≥ 6 earthquakes for the period 1988-present using modern broad-band seismograms. The main objective is to extend the database of energy-related parameters by combining 〈span〉mB〈/span〉 databases for recent and old events. The radiated energy 〈span〉ER_B〈/span〉 (in erg) computed from 〈span〉mB〈/span〉 using the Gutenberg & Richter relation log 〈span〉ER_B 〈/span〉=〈span〉 〈/span〉2.4〈span〉mB 〈/span〉+〈span〉 〈/span〉5.8 agrees very well with 〈span〉ER〈/span〉 estimated with modern techniques, especially for large deep earthquakes. Thus, 〈span〉ER_B〈/span〉 is useful as a proxy for 〈span〉ER〈/span〉 to investigate the global diversity of earthquake characteristics and physics over an extended period of time.〈/span〉

Description: 〈span〉〈div〉Summary〈/div〉We develop a new method for measuring ellipticity of Rayleigh waves from ambient noise records by degree-of-polarization (DOP) analysis. The new method, named DOP-E, shows a good capability to retrieve accurate ellipticity curves separated from incoherent noise. In order to validate the method we perform synthetic tests simulating noise in a 1D Earth model. We also perform measurements on real data from Antarctica and Northern Italy. Observed curves show a good fit with measurements from earthquake records and with theoretical ellipticity curves. The inversion of real data measurements for 〈span〉vS〈/span〉 structure shows a good agreement with previous models. In particular the shear-wave structure beneath Concordia station shows no evidence of a significant layer of liquid water at the base of the ice. The new method can be used to measure ellipticity at high frequency and therefore it will allow the imaging of near-surface structure, and possibly of temporal changes in subsurface properties. It promises to be useful to study near-surface processes in a wide range of geological settings, such as volcanoes, fault zones and glaciers.〈/span〉

Description: 〈span〉〈div〉Summary〈/div〉A gravity survey of the Al-Ain region within the United Arab Emirates (UAE) was designated to compare a regional gravity map with surface structural and geophysical studies of the area. The study area is located in the eastern part of the UAE near the Oman border. Gravity data were collected at 452 stations using a Scintrex CG-6 gravimeter. The gravity measurements define the subsurface density distributions in the study area. The Bouguer anomaly map shows values ranging from -21.5 mgal to 10.8 mgal. High Bouguer anomalies were observed in the eastern and western parts of the study area, in good agreement with the location of known thrust faults. Gravity gradient derivative methods (horizontal gradient, analytic signal, and tilt angle) were applied to the corrected gravity data to enhance visibility of subsurface geological features. The gradient derivative methods applied to the gravity data outline subsurface structures trending NNE-SSW, NE-SW and E-W. These trends can be correlated with regional tectonic phases and are in agreement with known surface structural trends. The power spectrum of the gravity data shows two distinct geological layers—a deeper layer from 2 km to 6 km depth representing Cretaceous formations, and a shallower layer from 0.5 km to 2 km depth representing Tertiary formations. The 3D gravity inversion results show subsurface density distribution characterized by high density contrast of geological bodies located in eastern and western parts of the study area, representing the Semail Ophiolites and Barzaman Formation, respectively. The Hafit Mountain is characterized by a low density contrast body at near surface (-500 kg/m〈sup〉3〈/sup〉) beneath which density increases with depth to + 50 kg/m〈sup〉3〈/sup〉 at 7 km (b.s.l.). This study is the first of its kind to image the regional subsurface and geological features of Al-Ain region with density distribution from ground surface to 15 km (b.s.l.).〈/span〉

Description: 〈span〉〈div〉Summary〈/div〉Deployment of cabled stations on the ocean bottom enhances the effectiveness of early warning systems for sub-oceanic earthquakes in terms of early signal detection and improved accuracy of hypocentre analyses. We investigated strong-motion data at cabled ocean-bottom stations deployed in the source area of the 2011 Tohoku earthquake (〈span〉M〈/span〉w 9.0) in northeastern Japan. We found a significantly large direct current component change in accelerations due to the roll and pitch motions of the cabled stations that are proportional to strong motions. We further found that the roll angle is affected more by the strong motions than the pitch angle owing to the cylindrical shape of the station housing in a long cable. The observed large roll motions caused anomalous displacement waveforms with significantly large amplitudes, which have contributed to the overestimation of magnitude in the two components perpendicular to the cable direction. Applying the scheme to three-component data for estimating real-time magnitude for early warnings in Japan resulted in high-magnitude estimates of M9.10–9.76 for the 2011 Tohoku earthquake, which are higher than the magnitude saturation and moment magnitude 〈span〉M〈/span〉w 9.0 of the earthquake. Such an overestimation was suppressed when the single-component data along the cable axis direction were used instead, with the resulting magnitude of M8.47–8.63 almost coinciding with the catalogue magnitude of 〈span〉M〈/span〉〈sub〉JMA〈/sub〉 8.4. For other event data that are not affected by the roll and pitch motions, no large systematic difference was found between the magnitudes estimated using the single- and three-component data. We therefore suggest that the use of single-component data could be an effective approach for stable magnitude estimation at cabled ocean-bottom stations for early warnings.〈/span〉

Description: 〈span〉〈div〉Summary〈/div〉An accurate traveltime approximation has a key role in the success of many seismic data processing, modeling, and inversion algorithms. In this paper, we introduce a new explicit six-parameter traveltime approximation as an extension to the known generalized moveout approximation. The parameters of the proposed approximation are estimated at two offsets. Zero-offset two-way time, normal moveout velocity, and effective anellipticity, are defined at zero offset; the remaining three parameters are defined from traveltime, ray parameter, and curvature at a reference offset. Each of these parameters have been previously employed in former traveltime approximations. The proposed method can be used to approximate full offset ray traced traveltimes, by two rays. Potential applications include approximation of migration Green's function, and common midpoint forward modeling. It also gives insight for traveltime behavior in anisotropic media. We provide numerical tests by different sources of nonhyperbolic moveout to show that the proposed method is a general improvement over the original generalized moveout approximation. It improves the accuracy, and makes the parameter selection more symmetrical, while adding a minimum computational burden.〈/span〉

Description: 〈span〉〈div〉Summary〈/div〉The recently evolving crustal deformation monitoring systems detect slow slip events, which are failures that are slower-paced than ordinary earthquakes, in the decollement zone of the shallow part of the plate subduction boundary zone. When the 2011 earthquake off the Pacific coast of Tohoku occurred, the decollement zone became the source of a massive tsunami. Deep-sea scientific drilling of the decollement zone has shown that the surrounding shear zone area is composed of a mixture of smectite and illite. These clay minerals differ in their frictional characteristics. This study proposes a mathematical model using a continuity equation to express the interactions between the solid and liquid phases and an elasto-plastic constitutive relation for the soil skeleton, considering the alteration of smectite to illite. In doing so, the decollement zone is viewed as a problem involving soil/water coupling in a mixture, where the soil skeleton is the solid and the pore water is the liquid. The dehydration of interlayer water in smectite, changes in frictional characteristics resulting from illitization, and changes in dilatancy characteristics are considered in this model. Shear simulations were also conducted, factoring in alteration. The simulation results revealed that the shear strength increased with alteration, and accordingly shear stress increased when the shear deformation is further applied. Since the shear strength increased with alteration, shear deformation was less likely to occur in the decollement zone, and it was considered that the sticking between the plate could be simulated.〈/span〉

Description: 〈span〉〈div〉Summary〈/div〉Crustal thickness (〈span〉H〈/span〉) and 〈span〉Vp〈/span〉/〈span〉Vs〈/span〉 ratio (κ) are two important parameters in characterizing a crustal structure. Their estimates from receiver function data often have large uncertainties due to the severe trade-off between the two parameters when the receiver functions only have the Moho 〈span〉P〈/span〉-to-〈span〉S〈/span〉 converted phase and its crustal multiples are weak or missing, which are common in tectonically active regions. Here, we extended the 〈span〉H〈/span〉-κ stacking method by combining receiver functions and 〈span〉PmP〈/span〉 travel times to determine 〈span〉H〈/span〉, κ, and 〈span〉Vp〈/span〉 simultaneously. We developed a deconvolution technique to measure 〈span〉PmP〈/span〉 travel times using the post-critical Moho 〈span〉P〈/span〉-wave reflection 〈span〉SsPmp〈/span〉 in the teleseismic 〈span〉S〈/span〉 wave. Synthetic tests show that the new methods can reliably extract the 〈span〉PmP〈/span〉 travel times and determine crustal thickness and velocities as long as the Moho 〈span〉Ps〈/span〉 is clear in receiver functions. We applied the methods to 63 broadband seismic stations in the western Hubei Province of China. The results confirm the rapid increase of crustal thickness revealed in a previous study from the low-elevation plains in the eastern part of the study area to the mountains in the west. The corresponding increase of crustal Poisson’s ratio from low-to-intermediate values to intermediate-to-high values suggests that the crust in the eastern part is more felsic, most likely due to delamination of a portion of the mafic lower crust of the Yangtze platform in the Late Mesozoic.〈/span〉

Description: 〈span〉〈div〉Abstract〈/div〉The purpose of this work is to generalize the frozen Gaussian approximation (FGA) theory to solve the 3-D elastic wave equation and use it as the forward modeling tool for seismic tomography with high-frequency data. FGA has been previously developed and verified as an efficient solver for high-frequency acoustic wave propagation (P-wave). The main contribution of this paper consists of three aspects: 1. We derive the FGA formulation for the 3-D elastic wave equation. Rather than standard ray-based methods (e.g. geometric optics and Gaussian beam method), the derivation requires to do asymptotic expansion in the week sense (integral form) so that one is able to perform integration by parts. Compared to the FGA theory for acoustic wave equation, the calculations in the derivation are much more technically involved due to the existence of both P- and S-waves, and the coupling of the polarized directions for SH- and SV-waves. In particular, we obtain the diabatic coupling terms for SH- and SV-waves, with the form closely connecting to the concept of Berry phase that is intensively studied in quantum mechanics and topology (Chern number). The accuracy and parallelizability of the FGA algorithm is illustrated by comparing to the spectral element method for 3-D elastic wave equation in homogeneous media; 2. We derive the interface conditions of FGA for 3-D elastic wave equation based on an Eulerian formulation and the Snell’s law. We verify these conditions by simulating high-frequency elastic wave propagation in a 1-D layered Earth model. In this example, we also show that it is natural to apply the FGA algorithm to geometries with non-Cartesian coordinates; 3. We apply the developed FGA algorithm for 3-D seismic wave-equation-based traveltime tomography and full waveform inversion, respectively.〈/span〉

Description: 〈span〉〈div〉Summary〈/div〉Accurate and efficient simulations of coseismic and postearthquake deformation are important for proper inferences of earthquake source parameters and subsurface structure. These simulations are often performed using a truncated halfspace model with approximate boundary conditions. The use of such boundary conditions introduces inaccuracies unless a sufficiently large model is used, which greatly increases the computational cost. To solve this problem, we develop a new approach by combining the spectral-element method with the mapped infinite-element method. In this approach, we still use a truncated model domain, but add a single outer layer of infinite elements. While the spectral elements capture the domain, the infinite elements capture the far-field boundary conditions. The additional computational cost due to the extra layer of infinite elements is insignificant. Numerical integration is performed via Gauss-Legendre-Lobatto and Gauss-Radau quadrature in the spectral and infinite elements, respectively. We implement an equivalent moment-density tensor approach and a split-node approach for the earthquake source, and discuss the advantages of each method. For postearthquake deformation, we implement a general Maxwell rheology using a second-order accurate and unconditionally stable recurrence algorithm. We benchmark our results with the Okada analytical solutions for coseismic deformation, and with the Savage & Prescott analytical solution and the PyLith finite-element code for postearthquake deformation.〈/span〉