The field of astronomy has arrived at a turning point in terms of size and complexity of both datasets and scientific collaboration. Commensurately, algorithms and statistical models have begun to adapt --- e.g., via the onset of artificial intelligence --- which itself presents new challenges and opportunities for growth. This white paper aims to offer guidance and ideas for how we can evolve our technical and collaborative frameworks to promote efficient algorithmic development and take advantage of opportunities for scientific discovery in the petabyte era. We discuss challenges for discovery in large and complex data sets; challenges and requirements for the next stage of development of statistical methodologies and algorithmic tool sets; how we might change our paradigms of collaboration and education; and the ethical implications of scientists' contributions to widely applicable algorithms and computational modeling. We start with six distinct recommendations that are supported by the commentary following them. This white paper is related to a larger corpus of effort that has taken place within and around the Petabytes to Science Workshops (https://petabytestoscience.github.io/).

The Photometric LSST Astronomical Time Series Classification Challenge (PLAsTiCC) is an open data challenge to classify simulated astronomical time-series data in preparation for observations from the Large Synoptic Survey Telescope (LSST), which will achieve first light in 2019 and commence its 10-year main survey in 2022. LSST will revolutionize our understanding of the changing sky, discovering and measuring millions of time-varying objects. In this challenge, we pose the question: how well can we classify objects in the sky that vary in brightness from simulated LSST time-series data, with all its challenges of non-representativity? In this note we explain the need for a data challenge to help classify such astronomical sources and describe the PLAsTiCC data set and Kaggle data challenge, noting that while the references are provided for context, they are not needed to participate in the challenge.

Galaxy mergers are key events in galaxy evolution, often causing massive starbursts and fueling active galactic nuclei (AGN). In these highly dynamic systems, it is not yet precisely known how much starbursts and AGN respectively contribute to the total luminosity, at what interaction stages they occur, and how long they persist. Here we estimate the fraction of the bolometric infrared (IR) luminosity that can be attributed to AGN by measuring and modeling the full ultraviolet to far-infrared spectral energy distributions (SEDs) in up to 33 broad bands for 24 merging galaxies with the Code for Investigating Galaxy Emission. In addition to a sample of 12 confirmed AGN in late-stage mergers, found in the $Infrared$ $Array$ $Satellite$ Revised Bright Galaxy Sample or Faint Source Catalog, our sample includes a comparison sample of 12 galaxy mergers from the $Spitzer$ Interacting Galaxies Survey, mostly early-stage. We perform identical SED modeling of simulated mergers to validate our methods, and we supplement the SED data with mid-IR spectra of diagnostic lines obtained with $Spitzer$ InfraRed Spectrograph. The estimated AGN contributions to the IR luminosities vary from system to system from 0% up to 91% but are significantly greater in the later-stage, more luminous mergers, consistent with what is known about galaxy evolution and AGN triggering.

Stars form in clustered environments, but how they form when the available resources are shared is still not well understood. A related question is whether the IMF is in fact universal across galactic environments, a galactic initial mass function (IGIMF), or whether it is an average of local IMFs. One of the long-standing problems in resolving this question and in the study of young clusters is observational: the emission from multiple sources is frequently seen as blended because at different wavelengths or with different telescopes the beam sizes are different. The confusion hinders our ability to fully characterize clustered star formation. Here we present a new method that uses a genetic algorithm and Bayesian inference to fit the blended SEDs and images of individual YSOs in confused clusters. We apply this method to the infrared photometry of a sample comprising 70 Spitzer-selected, low-mass ($M_{\rm{cl}}<100~\rm{M}_{\odot}$) young clusters in the galactic plane, and use the derived physical parameters to investigate the distributions of masses and evolutionary stages of clustered YSOs, and the implications of those distributions for studies of the IMF and the different models of star formation. We find that for low-mass clusters composed of class I and class II YSOs, there exists a non-trivial relationship between the total stellar mass of the cluster ($M_{\rm{cl}}$) and the mass of its most massive member ($m_{\rm{max}}$). The properties of the derived correlation are most compatible with the random sampling of a Kroupa IMF, with a fundamental high-mass limit of $150~\rm{M}_{\odot}$. Our results are also compatible with SPH models that predict a dynamical termination of the accretion in protostars, with massive stars undergoing this stopping at later times in their evolution.

We present the results of our investigation of the star-forming complexes W51 and W43, two of the brightest in the first Galactic quadrant. In order to determine the young stellar object (YSO) populations in W51 and W43 we used color-magnitude relations based on Spitzer mid-infrared and 2MASS/UKIDSS near-infrared data. We identified 302 Class I YSOs and 1178 Class II/transition disk candidates in W51, and 917 Class I YSOs and 5187 Class II/transition disk candidates in W43. We also identified tens of groups of YSOs in both regions using the Minimal Spanning Tree (MST) method. We found similar cluster densities in both regions even though Spitzer was not able to probe the densest part of W43. By using the Class II/I ratios, we traced the relative ages within the regions and based on the morphology of the clusters we argue that several sites of star formation are independent of one another in terms of their ages and physical conditions. We used spectral energy distribution (SED)-fitting to identify the massive YSO (MYSO) candidates since they play a vital role in the star formation process and then examined them to see if they are related to any massive star formation tracers such as UCH II regions, masers and dense fragments. We identified 17 MYSO candidates in W51, and 14 in W43, respectively and found that groups of YSOs hosting MYSO candidates are positionally associated with H II regions in W51, though we do not see any MYSO candidates associated with previously identified massive dense fragments in W43.

(abridged) The majority of star-forming galaxies follow a simple empirical correlation in the star formation rate (SFR) versus stellar mass ($M_*$) plane, usually referred to as the star formation Main Sequence (MS). Here we combine a set of hydro-dynamical simulations of interacting galactic disks with state-of-the-art radiative transfer codes to analyze how the evolution of mergers is reflected upon the properties of the MS. We present \textsc{Chiburst}, a Markov Chain Monte Carlo (MCMC) Spectral Energy Distribution (SED) code that fits the multi-wavelength, broad-band photometry of galaxies and derives stellar masses, star formation rates, and geometrical properties of the dust distribution. We apply this tool to the SEDs of simulated mergers and compare the derived results with the reference output from the simulations. Our results indicate that changes in the SEDs of mergers as they approach coalescence and depart from the MS are related to an evolution of dust geometry in scales larger than a few hundred parsecs. This is reflected in a correlation between the specific star formation rate (sSFR), and the compactness parameter $\mathcal{C}$, that parametrizes this geometry and hence the evolution of dust temperature ($T_{\rm{dust}}$) with time. As mergers approach coalescence, they depart from the MS and increase their compactness, which implies that moderate outliers of the MS are consistent with late-type mergers. By further applying our method to real observations of Luminous Infrared Galaxies (LIRGs), we show that the merger scenario is unable to explain these extreme outliers of the MS. Only by significantly increasing the gas fraction in the simulations are we able to reproduce the SEDs of LIRGs.

Spatially resolved kinematics have been used to determine the dynamical status of star-forming galaxies with ambiguous morphologies, and constrain the importance of galaxy interactions during the assembly of galaxies. However, measuring the importance of interactions or galaxy merger rates requires knowledge of the systematics in kinematic diagnostics and the visible time with merger indicators. We analyze the dynamics of star-forming gas in a set of binary merger hydrodynamic simulations with stellar mass ratios of 1:1 and 1:4. We find that the evolution of kinematic asymmetries traced by star-forming gas mirrors morphological asymmetries derived from mock optical images, in which both merger indicators show the largest deviation from isolated disks during strong interaction phases. Based on a series of simulations with various initial disk orientations, orbital parameters, gas fractions, and mass ratios, we find that the merger signatures are visible for ~0.2-0.4 Gyr with kinematic merger indicators but can be approximately twice as long for equal-mass mergers of massive gas-rich disk galaxies designed to be analogs of z~2-3 submillimeter galaxies. Merger signatures are most apparent after the second passage and before the black holes coalescence, but in some cases they persist up to several hundred Myr after coalescence. About 20-60% of the simulated galaxies are not identified as mergers during the strong interaction phase, implying that galaxies undergoing violent merging process do not necessarily exhibit highly asymmetric kinematics in their star-forming gas. The lack of identifiable merger signatures in this population can lead to an underestimation of merger abundances in star-forming galaxies, and including them in samples of star-forming disks may bias the measurements of disk properties such as intrinsic velocity dispersion.

The total infrared (IR) luminosity is very useful for estimating the star formation rate (SFR) of galaxies, but converting the IR luminosity into an SFR relies on assumptions that do not hold for all galaxies. We test the effectiveness of the IR luminosity as an SFR indicator by applying it to synthetic spectral energy distributions generated from three-dimensional hydrodynamical simulations of isolated disc galaxies and galaxy mergers. In general, the SFR inferred from the IR luminosity agrees well with the true instantaneous SFR of the simulated galaxies. However, for the major mergers in which a strong starburst is induced, the SFR inferred from the IR luminosity can overestimate the instantaneous SFR during the post-starburst phase by greater than two orders of magnitude. Even though the instantaneous SFR decreases rapidly after the starburst, the stars that were formed in the starburst remain dust-obscured and thus produce significant IR luminosity. Consequently, use of the IR luminosity as an SFR indicator may cause one to conclude that post-starburst galaxies are still star-forming, whereas in reality, star formation was recently quenched.

A significant fraction ($\sim 30$ %) of the gamma-ray sources listed in the second $\textit{Fermi}$ LAT (2FGL) catalog is still of unknown origin, being not yet associated with counterparts at lower energies. Using the available information at lower energies and optical spectroscopy on the selected counterparts of these gamma-ray objects we can pinpoint their exact nature. Here we present a pilot project pointing to assess the effectiveness of the several classification methods developed to select gamma-ray blazar candidates. To this end, we report optical spectroscopic observations of a sample of 5 gamma-ray blazar candidates selected on the basis of their infrared WISE colors or of their low-frequency radio properties. Blazars come in two main classes: BL Lacs and FSRQs, showing similar optical spectra except for the stronger emission lines of the latter. For three of our sources the almost featureless optical spectra obtained confirm their BL Lac nature, while for the source WISEJ022051.24+250927.6 we observe emission lines with equivalent width $EW\sim 31$ $\AA$, identifying it as a FSRQ with $z = 0.48$. The source WISEJ064459.38+603131.7, although not featuring a clear radio counterpart, shows a blazar-like spectrum with weak emission lines with $EW \sim 7$ $\AA$, yielding a redshift estimate of $z=0.36$. In addition we report optical spectroscopic observations of 4 WISE sources associated with known gamma-ray blazars without a firm classification or redshift estimate. For all of these latter sources we confirm a BL Lac classification, with a tentative redshift estimate for the source WISEJ100800.81+062121.2 of $z = 0.65$.

NGC 604 is the second most massive H II region in the Local Group, thus an important laboratory for massive star formation. Using a combination of observational and analytical tools that include Spitzer spectroscopy, Herschel photometry, Chandra imaging, and Bayesian Spectral Energy Distribution fitting, we investigate the physical conditions in NGC 604, and quantify the amount of massive star formation currently taking place. We derive an average age of 4 +/- 1 Myr and a total stellar mass of 1.6 (+1.6)(-1.0) x 10^5 M_sun for the entire region, in agreement with previous optical studies. Across the region we find an effect of the X-ray field on both the abundance of aromatic molecules and the [Si II] emission. Within NGC 604 we identify several individual bright infrared sources with diameters of about 15 pc and luminosity weighted masses between 10^3 M_sun and 10^4 M_sun. Their spectral properties indicate that some of these sources are embedded clusters in process of formation, which together account for ~8% of the total stellar mass in the NGC 604 system. The variations of the radiation field strength across NGC 604 are consistent with a sequential star formation scenario, with at least two bursts in the last few million years. Our results indicate that massive star formation in NGC 604 is still ongoing, likely triggered by the earlier bursts.

To understand and interpret the observed Spectral Energy Distributions (SEDs) of starbursts, theoretical or semi-empirical SED models are necessary. Yet, while they are well-founded in theory, independent verification and calibration of these models, including the exploration of possible degeneracies between their parameters, are rarely made. As a consequence, a robust fitting method that leads to unique and reproducible results has been lacking. Here we introduce a novel approach based on Bayesian analysis to fit the Spitzer-IRS spectra of starbursts using the SED models proposed by Groves et al. (2008). We demonstrate its capabilities and verify the agreement between the derived best fit parameters and actual physical conditions by modelling the nearby, well-studied, giant HII region 30 Dor in the LMC. The derived physical parameters, such as cluster mass, cluster age, ISM pressure and covering fraction of photodissociation regions, are representative of the 30 Dor region. The inclusion of the emission lines in the modelling is crucial to break degeneracies. We investigate the limitations and uncertainties by modelling sub-regions, which are dominated by single components, within 30 Dor. A remarkable result for 30 Doradus in particular is a considerable contribution to its mid-infrared spectrum from hot ({\simeq} 300K) dust. The demonstrated success of our approach will allow us to derive the physical conditions in more distant, spatially unresolved starbursts.

We present a model that describes stellar infrared excesses due to heating of the interstellar (IS) dust by a hot star passing through a diffuse IS cloud. This model is applied to six lambda Bootis stars with infrared excesses. Plausible values for the IS medium (ISM) density and relative velocity between the cloud and the star yield fits to the excess emission. This result is consistent with the diffusion/accretion hypothesis that lambda Bootis stars (A- to F-type stars with large underabundances of Fe-peak elements) owe their characteristics to interactions with the ISM. This proposal invokes radiation pressure from the star to repel the IS dust and excavate a paraboloidal dust cavity in the IS cloud, while the metal-poor gas is accreted onto the stellar photosphere. However, the measurements of the infrared excesses can also be fit by planetary debris disk models. A more detailed consideration of the conditions to produce lambda Bootis characteristics indicates that the majority of infrared-excess stars within the Local Bubble probably have debris disks. Nevertheless, more distant stars may often have excesses due to heating of interstellar material such as in our model.

We have discovered a bow shock shaped mid-infrared excess region in front of delta Velorum using 24 micron observations obtained with the Multiband Imaging Photometer for Spitzer (MIPS). The excess has been classified as a debris disk from previous infrared observations. Although the bow shock morphology was only detected in the 24 micron observations, its excess was also resolved at 70 micron. We show that the stellar heating of an ambient interstellar medium (ISM) cloud can produce the measured flux. Since delta Velorum was classified as a debris disk star previously, our discovery may call into question the same classification of other stars. We model the interaction of the star and ISM, producing images that show the same geometry and surface brightness as is observed. The modeled ISM is 15 times overdense relative to the average Local Bubble value, which is surprising considering the close proximity (24 pc) of delta Velorum. The abundance anomalies of lambda Bootis stars have been previously explained as arising from the same type of interaction of stars with the ISM. Low resolution optical spectra of delta Velorum show that it does not belong to this stellar class. The star therefore is an interesting testbed for the ISM accretion theory of the lambda Bootis phenomenon.