Quantifying burstiness, and another brief look at SDSS J095343.89-000524.7

One simple way to quantify the burstiness of star formation is just to estimate the average star formation rate over large time intervals divided by the average SFR over cosmic time. Of particular interest is the time interval between ~100 Myr and ~1 Gyr since this is roughly the time interval that a post-starburst galaxy is recognizable as such.

Partly because it happens to still be in my active workspace and partly because it’s really interesting I’m going to take another look at SDSS J095343.89-000524.7 (MaNGA mangaid 1-897).  This was in the post-starburst ancillary sample, selected from the catalog by Pattarakijwanich et al.

This image from the Subaru HSC-SSP survey1retrieved as a screenshot from the Legacy Survey sky browser. is much deeper than SDSS imaging and clearly shows extended tidal tails and debris, suggesting that these galaxies have been interacting for some time.

SDSS J095343.89-000524.7 (observed as mangaid 1-897). Image screenshot from Subaru HSC survey.

Moving on to various properties derived from the MaNGA spectroscopy and my SFH models with, still, EMILES based SSP models. First here are maps of stellar mass density and 100 Myr averaged star formation rate density. Note that I rebinned the spectra from two posts ago to try to capture more of the tidal tails while excluding the truly blank regions of sky. There are two clear peaks in the stellar mass density separated by a projected distance of about 11 kpc. The central stellar mass densities are nearly the same at about 108.95 M☉/kpc2 . Interestingly enough the bright white peak in surface brightness appears not to coincide with the western peak in stellar mass density, but is offset by a small amount to the north.

Note also that the highest recent star formation is offset to the north of the apparent western nucleus. I’ll look at that in more detail below.

MaNGA plateifu 10843-9101 (mangaid 1-897). Maps of stellar mass density and star formation rate density.

The ionized gas properties are rather different in the two galaxies. Below are BPT classifications using, as usual for me, just the [O III]/Hβ vs. [N II]/Hα diagnostics and Kauffmann’s classification scheme. Emission line fluxes are generally stronger in the eastern galaxy with mostly star forming line ratios. Note two spectra with “composite” line ratios are near the eastern nucleus and might therefore actually be due to a mix of stellar and AGN ionization.

MaNGA plateifu 10843-9101 (mangaid 1-897). BPT classifications from [O III]/Hβ vs. [N II]/Hα diagnostics

I calculate a few “strong line” gas metallicity estimates from standard literature sources. The one that seems to produce the most consistent estimates is the calibration of Dopita et al. (2016) based on the ratios of [N II 6548]/[S II 6717, 6731] and [N II]/Hα. The eastern galaxy shows a fairly smooth radial gradient while the west is considerably metal enriched in the region with the strongest starburst. The highest metallicity is right at the center of the IFU at the position of the bright white source.

MaNGA mangaid 1-897 (plateifu 10843-9101). Gas phase metallicity 12 + log(O/H) from strong line calibration of Dopita et al. (2016).

Let’s return to the idea I had at the top of the post to look at star formation rates in broad time intervals relative to the mean star formation rate over cosmic time. For this exploratory exercise I used just 4 bins with upper age limits of 0.1, 1.25, 2.25, and (nominally) 14 Gyr. There seems no point being too fastidious about calculating the bin widths: I just used the difference in nominal ages between the endpoints. I did take into account the lookback time to the galaxy, which for this one is about 1 Gyr (z = 0.083), so the final bin has a calculated width of 10.5Gyr. I chose to make the 3rd, intermediate age bin a rather short 1 Gyr wide to look for aging starbursts that might be missed using the typical selection criterion of strong Balmer absorption. In this case there’s no evidence of that: both galaxies seem to have had uneventful histories up until ~1 Gyr ago.

The top row of the plot below is the most interesting: there appear to have been two major bursts of recent star formation, both highly localized to the central region of the western galaxy. If the model estimate of the location of the peak stellar mass density is correct the fiber with the largest star formation excess in the 100 Myr – 1.25Gyr interval is offset just to the north and coincident with the IFU center. The more recent burst is also offset from the older one. There is a hint of recent accelerated star formation over most of both galaxies.

MaNGA plateifu 10843-9101 (mangaid 1-897). Maps of relative average SFR over the designated time intervals.

For the rest of this post I plot model fits to the spectra and star formation histories for the fibers surrounding the two nuclei. These are ordered approximately from north to south and west to east. For reference the IFU center is at (ra, dec) = (148.43291, -0.09018). The model has the peak stellar mass density in the western system at (ra, dec) = (148.4328, -0.09062). The eastern galaxy’s nucleus is at (ra, dec) = (148.4349, -0.09064).

Note below that the plots have different vertical scales. The horizontal scales are the same for both spectra and star formation histories, but at least one SFH plot is slightly misaligned.

Central region – western galaxy

Central region – eastern galaxy

In an earller post I mentioned a MaNGA related paper by Cheng et al. who found nearly 500 systems with post-starburst characteristics that fell in 3 broad categories: centrally concentrated PSB regions, ring-like, and irregularly located. Clearly any galaxy that was selected based on SDSS spectroscopy that’s not a false positive will have a central PSB region, although that of course doesn’t preclude extended post-starburst conditions. This particular galaxy appears to have a remarkably compact post-starburst region.

When time permits again I plan to look at the remaining 40 galaxies in this sample. Unfortunately the larger sample of Cheng et al. appears to have no published catalog.

UGC 10200 and MCG +07-33-040

The Hubble Space Telescope “gap filler” program “Gems of the Galaxy Zoos” (proposal ID 15445, PI William Keel) had several prospective targets that I played a small role in selecting, and this recent HST observation was one of them. The actual target was the small disturbed galaxy at top left, which I will refer to as MCG +07-33-040. I don’t know if it was fortuitous that the larger and brighter UGC 10200 was also imaged in the same ACS field, but these are clearly interacting or at least have in the recent past, as is the small system in the upper right, which is identified as a blue compact galaxy with redshift z=0.00556 in Pustilnik et al. (1999). I’m going to focus on the top left galaxy in this post.

Galaxies UGC 10200 (lower right) and MCG +07-33-040 (upper left). HST/ACS, F475W filter. Proposal ID 15445, PI Keel.

What interested me wasn’t the galaxy image so much as its SDSS spectrum, which has three interesting characteristics:

SDSS spectrum of central part of MCG +07-33-040

First, this is a classic post starburst galaxy spectrum with extremely strong Balmer absorption lines1My code measures the Lick index HδA as an exceptionally strong 8.06 ± 0.41 Å. and no obvious evidence of emission. In fact, although this designation isn’t used much anymore, it’s actually a classic “A+K” spectrum which reverses the usual “K+A” terminology to indicate the light is dominated by early type (i.e. young) stars. Second and third, the spectrum was misclassified as coming from a white dwarf star, and the redshift was erroneously estimated as around 0.004 which was the maximum allowed for stars in the SDSS data reduction pipeline. Using a variation of the code that I use to measure redshift offsets I get a robust value of z = 0.006682 ± 9E-06

Template fit to SDSS spectrum of MCG +07-33-040

This is almost exactly the same redshift as its nearby companion UGC 10200 (also in the HST image above), which has a securely determined z = 0.00664

SDSS spectrum of central region of UGC 10200

For the rest of this post I’m going to assume the Hubble flow redshift is the measured one, which with my adopted cosmological parameters2which for the record are H0 = 70 km/sec/Mpc, Ωm = 0.27, Ωλ = 0.73. make the luminosity distance 28.8 Mpc, the distance modulus m-M = 32.3 mag, and the angular scale 138 pc/” or about 7 pc per ACS pixel. The projected distance between the centers of the two bright galaxies in the HST image is about 96″ or 13.2 kpc.

I spent some time last weekend downloading and starting to learn the software Aperture Photometry Tool (APT), which is interactive software for manually performing aperture photometry. Zooming in on the center of the presumed post starburst galaxy I located the reported position of the SDSS fiber as closely as I could. In the screenshot below the aperture radius was set to 30 pixels, the same size as the SDSS spectroscopic fibers. I measured the F475W AB magnitude to be 17.90 ± 0.013 without sky subtraction, which is close enough to the SDSS g band fiberMag estimate of 18.05. The SDSS g band Petrosian magnitude estimate is 15.16, so the fiber contains about 7% of the total galaxy light.

Central region of MCG +07-33-040 with position and size of SDSS fiber overlaid. Screenshot from APT

A striking feature of the HST image is there are many point-like symmetrical objects embedded within the otherwise nearly featureless diffuse light of the galaxy. Within the SDSS fiber footprint I count about 8-10 of these (the range being due to some uncertainty about what to call point-like and symmetrical). In order to get a handle on their contribution to the spectrum I did aperture photometry on them using a 3 pixel radius aperture with median sky subtraction from a 5 to 8 pixel radius annulus. The apparent magnitudes of the 5 brightest objects range from about 22.6 to 23.1. The summed luminosity of those 5 amounts to only 3.5% of the total light in the fiber, so the spectrum is mostly telling us something about the diffuse starlight. Even if one or more of those objects are foreground stars they can’t be a significant source of contamination. Clicking around the blank regions of the HST field I found fewer than one star per SDSS fiber size region, so it’s likely there are few if any foreground stars within the visible extent of the galaxy.

There is plenty of observational and theoretical evidence that massive star clusters are formed in mergers and close encounters of galaxies. As a coincidental example the merger remnant NGC 3921 — which was one of the 4 galaxies in my last post — has over 100 young globular clusters located both in the main body and southern tidal tail (Schweizer et al. 1996; Knierman et al. 2003). The brightest source in this galaxy (near the southern edge of the visible fuzz) has an apparent magnitude of m ≈ +21.7, which for the adopted distance modulus is M ≈ -10.6. With a solar g band absolute magnitude of 5.11 this corresponds to L ≈ 1.9×106 L . The 5 brightest objects within the fiber have absolute magnitudes between about -9.7 and -9.2. These would be quite luminous for galactic globular clusters, but they’re likely to be fairly young and would fade by at least a few magnitudes as they age.

I haven’t tried a more sophisticated analysis of these objects’ sizes, but using the APT radial profile tool the presumed clusters look little different from nearby foreground stars and all that I’ve examined have FWHM diameters around 2-2.5 pixels. A strict upper limit to their sizes is therefore around 14 pc.

Someday I may undertake a complete census and luminosity function of the cluster system in this galaxy, and perhaps also look at the neighboring starburst galaxy UGC 10200. These systems by the way are cataloged as an interacting dwarf galaxy pair by Paudel et al. (2018) with a total stellar mass of log(M*) = 9.5 and a 3:1 mass ratio, which makes the estimated stellar mass of this galaxy just under 109 M. The system is very gas rich, with a neutral hydrogen mass estimated (by the same source) of 109.69 M. There are actually at least two published HI maps of this system. The one below, from Thomas et al. (2004) shows that atomic hydrogen extends over essentially the entire extent of the Hubble image above, including the target galaxy.

VLA map of HI gas in UGC 10200 system

Next I turn to star formation history models for the post starburst spectrum at the top of the post. This uses the same Stan model code as my MaNGA investigations with some minor pre- and post-processing adjustments. I ran two separate models. One used a metal poor subset of the EMILES SSP libraries with Z ∈ {[-2.27], [-1.26], [-0.25]} with, as usual, Kroupa IMF and BaSTI isochrones. I did not attempt to append younger models, so the youngest age is 30Myr. For completeness I also ran a model with my usual EMILES subset + PYPOPSTAR models and Z ∈ {[-0.66], [-0.25], [+0.06], [+0.40]}. First, here is the modeled star formation history with the metal poor subset. I’ve again used a logarithmic time scale and linear star formation rate scale.

Model star formation history of central region of MCG +07-33-040 using metal poor subset of EMILES SSP library

Next is the metal rich subset:

Model star formation history of central region of MCG +07-33-040 using metal rich subset of EMILES+pypopstar SSP library

Both model runs show a fairly steep ramp up in star formation beginning at about 600Myr lookback time and a steep decline around 50Myr ago. The lingering star formation in the metal rich model might be a manifestation of the infamous “age metallicity degeneracy” since Balmer Hα emission is too low to support this level of star formation. Comparing the mass growth histories exposes a more subtle effect: the metal poor models have a consistently higher mass fraction at any given epoch. Also, the period of accelerated star formation involved a slightly smaller fraction of the present day stellar mass.

Mass growth histories of MCG +07-33-040 using metal poor and metal rich subsets of EMILES SSP library

Both models fit the data well. In terms of mean log-likelihood the metal poor model outperformed the metal rich, but only by about 0.4%. The proper Bayesian way to compare models is through the “evidence,” which is hard to estimate accurately. I suspect the metal poor model would be at least slightly flavored because it has fewer parameters than the metal rich one.

Posterior predictive fit to SDSS spectrum of MCG +07-33-040

The duration of accelerated star formation (about which both models agree) is a little surprising in light of simulations that usually show a fairly short SF burst in the first passage in mergers. But, simulations have only explored a small range of the potential parameter range. Studies of low mass galaxies with extended, massive HI haloes might be of interest.

One more sanity check. Suppose the closest approach between our target and UGC 10200 was 60Myr ago, allowing another 10Myr before (presumably) supernova feedback quenched star formation. Assuming the relative motion is transverse to our line of sight traveling 13.2 kpc in 60Myr implies an average separation speed of ≈215 km/sec. This is a perfectly reasonable value for a galaxy pair or loose group.

Finally for this spectrum, here is a quick look at emission line fluxes. Even though visually not at all obvious several lines were detected at marginal (>2σ) to high (>5σ) confidence. A couple of surprises are the [O I] 6300Å line, which is often only marginally detected even in star forming systems, is a firm (3σ) detection and stronger than the usually more prominent [O III] doublet. Also, the [S II] 6717-6730 doublet is a firm detection while the [N II] doublet is not. What this means is unclear to me. Most of the “strong emission line” metallicity indicators that I have formulae for include [N II] (or [O II] which are out of the wavelength range of these spectra), so it isn’t really possible to make a gas metallicity estimate. It’s a safe guess it’s subsolar though.

line[Ne III] 3869[Ne III] 3970[O III] 4959[O III] 5007[O I] 6300[O I] 6363[N II] 6548[N II] 6584[S II] 6717[SII] 6730
mean17.12.31.51.61.92.17.92.44.98.22.82.939.12.514.414.2
s.d.6.32.01.41.41.61.83.12.02.92.81.92.02.61.82.82.8
ratio2.71.11.11.11.21.22.61.21.73.01.51.515.21.45.25.2
Flux measurements for tracked emission lines in spectrum of MCG +07-33-040. Units are 10-17 erg/sec/cm2

There are at least two questions about this galaxy that it would be nice to have answers for. First, since the SDSS fiber only includes about 7% of the luminosity and a similar fraction of the stellar mass we really don’t know if it is recently quenched globally or just where SDSS happened to target. My guess from this HST image is that it is globally quenched because its companion UGC 10200 shows clear evidence of dust lanes and extended star forming regions even in this monochromatic image, while the diffuse light in this galaxy looks relatively featureless. A definitive answer would require IFU spectroscopy though.

A second question is whether the star cluster system is truly young or primordial (or both). This would require color measurements from a return visit by HST using at least one more filter in the red. Estimating a luminosity function is feasible with the existing data, although it would have rather shallow coverage. From my casual clicking around the image it appears to be possible to reach magnitudes a little larger than +24 with reasonable precision.

When this topic first came up on the old Galaxy Zoo talk I thought these might comprise a new and overlooked category of galaxies. In fact though all of the examples I investigated belonged to cataloged galaxies and most of the spectra were of small regions in much larger nearby galaxies. A few galaxies that were in the original Virgo Cluster Catalog and excluded from the EVCC because of lack of redshift confirmation should be added back. There were probably only a few like this one with large errors in redshift estimates and high signal to noise spectra. I haven’t spent enough time with the literature to know if rapidly quenched dwarf galaxies are especially interesting. Maybe they are.

Journal notes: Haines et al. (2015), “Testing the modern merger hypothesis…”

While browsing through the ADS listing of papers that cite Schawinski’s paper that I’ve been discussing for a while I came across this one by Haines et al. with the full title “Testing the modern merger hypothesis via the assembly of massive blue elliptical galaxies in the local Universe”. Besides being on the same theme of searching for post-starburst or “transitional” galaxies in the local universe that I’ve been pursuing for some time the paper was interesting because it made use of IFU based spectroscopic data that predates MaNGA. As it happens 4 of the 12 galaxies have observations in the final MaNGA release, providing an excellent opportunity to compare results from completely independent data sets.

The “modern merger hypothesis” that the authors tested relates to a topic I’ve discussed before, which is that N-body simulations show that strong, centrally concentrated starbursts are a possible outcome of major gas rich galaxy mergers around the time of coalescence. If some feedback process (an AGN or supernovae) rapidly quenches star formation there will ensue a period of time when the galaxy will be recognizable as post-starburst.

In a series of long and rather difficult (and influential judging by the number of citations) Hopkins and collaborators (2006, 2008a, 2008b) have made a case that major gas rich mergers with accompanying starbursts are in fact the major pathway to the formation of modern elliptical galaxies. They claim that their merger hypothesis accounts for a variety of phenomena, including the growth and evolution of supermassive black holes and quasars.

The specific aspect of the merger hypothesis this study tried to address was the prevalence of strong centrally concentrated starbursts in a sample of ellipticals in the process of forming as evidenced by visible disturbances consistent with recent mergers. The main tool they used was a suite of simple star formation history models with exponentially decaying star formation rate with single (also exponentially decaying) bursts on top of varying ages and decay time scales. They used these to predict just two quantities: Balmer absorption line strength measured by the average of the Lick HδA and HγA indexes, and the 4000Å break strength index Dn4000. For reference here is a screen grab of their model trajectories:

Predected trajectories in the Hδ – Dn4000 plane per Haines et al. (2015). Clipped from the electronic journal paper.

This is a pretty standard calculation variations of which have been performed for decades, and this graph looks much like others I have seen in the literature. A fairly basic problem with it though is that position in the Balmer – D4000 plane doesn’t uniquely constrain even the recent stellar evolution. In astronomers’ parlance there is a “degeneracy”1the term refers to a situation in which multiple combinations of some parameters of interest produce effectively equivalent values of some observable(s), or of course the converse. The best known example is the “age-metallicity degeneracy,” which refers to the fact that an old metal poor population looks like a younger metal rich one in several respects such as broad band colors. between burst strength (if any) and burst age. This is a well known problem with the Balmer line strength index that was already recognized by Worthey and Ottaviani (1997), who developed these indexes. Adding a second index in the form of the 4000Å break strength doesn’t break the degeneracy: there are regions of the plane where bursting and non-bursting populations overlap, as can be seen clearly in the graphic above. This is actually a problem for any attempt to identify post-starburst galaxies. After correcting for emission most ordinary starforming galaxies have strong Balmer absorption lines, so using that index alone will certainly produce many false positives. On the other hand selection criteria like those used by Goto and many others before and after — selecting for both strong Balmer absorption and weak emission — will capture only a small interval in post-starburst galaxies’ life cycles.

hd_d4000_bigsample
Hδ line strength vs. 4000Å break index for a large (~380K) sample of SDSS galaxy spectra. Measurements from the MPA-JHU analysis pipeline downloaded from SDSS Skyserver

Let’s get to results. Some basic details of the sample are in the table below. Morphological classifications are from McIntosh et al. (2014) as given in this paper. The abbreviations are SPM: spherical post merger; pE: peculiar Elliptical. The two marked pE/SPM didn’t have a strong consensus among several professional classifiers. I list them in order of my own visual impression of degree of disturbance. I also list redshifts taken from the MaNGA catalog and Petrosian colors.

NED nameNYU IDmangaidplateifuMorphzu-rg-i
NGC 39215410441-61744510510-6103SPM0.0191.970.86
MRK 3857194861-6049708940-6102pE/SPM0.0281.430.63
MRK 3661009171-6033097993-1902pE/SPM0.0271.590.79
NGC 1149223181-371558154-6103pE0.0292.291.11
Columns: (1) Common catalog designation (NED name). (2) NYU VAC ID. (3) MaNGA mangaid. (4) MaNGA plateifu. (5) Morphology (see text). (6) redshift from MaNGA DRP catalog. (7-8) Petrosian u-r and g-i colors from NYU VAC via the MaNGA DRP catalog.

The main prediction of the merger with accompanying centrally concentrated starburst hypothesis the paper tests is that the Balmer absorption index should be large and have a negative gradient with radius while the 4000Å break strength should be low with a positive gradient. The authors concluded that only one member of their sample — nyu541044 — clearly falls in the post-starburst region (marked as region 4 in the graph above) of the <Hδ, Hγ> – Dn4000 plane. The two pE/PM galaxies, both of which are in my sample, lie in the starforming region 1. They inferred from this that these galaxies are undergoing at most a weak burst. I’m going to mildly disagree with that conclusion.

Screenshot from 2022-07-07 15-23-36
Measured values for the specified indexes from Haines et al. (2015). Clipped from the electronic journal paper.

I have calculated the pseudo Lick index HδA and Dn4000 as part of my analysis “pipeline” since I started this hobby. I actually make these measurements in the initial maximum likelihood fitting step since they don’t depend on modeling except for small (usually) emission corrections. I don’t calculate an Hγ index, but its theoretical behavior is similar to Hδ. I’m trying here just to verify the approximate magnitude and radial trends of the chosen indexes. The two IFUs used in the Haines study had larger spatial coverage than these MaNGA observations (but much smaller wavelength coverage, which will become important). Instead of their strategy of binning in annuli I used my usual Voronoi binning strategy with a minimum target S/N. There were some oddities in the NYU estimates of effective radii so I chose to use distances from the IFU center in kpc for these plots. The distances assigned to the multiply binned spectra are the same as Cappelari’s published code produces; for single fiber spectra it’s just the position of the fiber center.

My measurements agree reasonably well with those of Haines et al. All three of the most disturbed galaxies have central Hδ indexes > 5Å with NGC 3921 (plateifu 10510-6103, nyu541044) having a larger central value and steeper gradient in the inner few kpc than the two pE/SPM galaxies. The fourth galaxy shows no obvious trend in either index with radius2The next several plots show trend lines for each galaxy computed by fitting simple loess curves to the data using the default parameters in ggplot2. These, and especially the confidence bands included in the plots, should not be taken seriously!. The central values where the S/N is highest are in good agreement.

Lets turn to the results of star formation history models, which I ran on all 4 data sets. First, here are 100Myr averaged star formation rate density and specific star formation rate versus distance:

Star formation rate density vs. distance from IFU center (kpc) for 4 disturbed early type galaxies.
Specific star formation rate density vs. distance from IFU center (kpc) for 4 disturbed early type galaxies.

Three of these galaxies are clearly experiencing centrally concentrated episodes of star formation, and two are at or near starburst levels in specific star formation rate near their centers. As seen below two of these straddle my estimate of the “spatially resolved star forming main sequence” while the one presumed post-starburst galaxy reaches it in the central region.

mstar_sfr_4spm
Star formation rate density versus stellar mass density for 4 disturbed early type galaxies

As I’ve shown several times before there’s a reasonably tight linear relationship between modeled star formation rate and Hα luminosity density. The plot shows Hα luminosity density corrected for modeled stellar redenning, which certainly underestimates attenuation in emission regions. The modeled star formation rates are consistently above the Kennicut relation shown as the straight line as I’ve seen in every sample I’ve looked at.

Star formation rate density vs. Hα luminosity density for 4 disturbed early type galaxies

Finally, lets take a look at detailed star formation histories. Instead of my usual practice of plotting them all in a grid here I just display 2 binned star formation histories. One comprises the innermost 7 bins, which since the fibers are arranged in a hexagonal grid should form a regular hexagon around the IFU center. These range in “radius” from about 0.75 to 1.1 kpc in these four galaxies. The second is for an “annulus” in approximately the outer kpc of each IFU. The extent of the IFU footprints ranges from 3.1 to 5.9 kpc. I calculate these by summing the contributions in each SFH model contributing to the bins, not by running new models for binned spectra. Since the dithered fiber positions overlaps this overestimates the total mass in each bin, but I care about the shape and timing of events rather than the absolute values of star formation rate estimates.

The next 4 plots display the results. Lookback time is logarithmically scaled with the same range and ticks for each SFH. Vertical scales are linear and differ for each graph. The graphs are in the same order as the basic information table above. As I’ve written before these models “want” to have smoothly varying mass per time bin which has the unfortunate effect of producing jumps in the apparent SFR when the bin widths change. In the BaSTI isochrone based SSP models these occur at 100 Myr, 1 Gyr, and 4 Gyr and can sometimes be quite prominent.

With caveats out of the way the one clear post-starburst in the sample had (per the model) a powerful and short starburst at ≈300 Myr lookback time, with a small amount continuing to the present (this can’t be seen at the scale of the graph, but ongoing star formation is ~1 M/yr). The total mass contribution from the burst and subsequent star formation is around 15%.

The two apparent ongoing starbursts have later bursts of star formation that are slightly weaker in terms of total mass contribution and peak star formation rate, but still quite significant. All three of the starburst/post-starburst galaxies appear to have had two major waves of late time (last ~2 Gyr or less) star formation. As I’ve written before in merger simulations the progenitors usually complete a few orbits before coalescence, with some enhanced star formation around each perigalactic passage. I hesitate to take these models that literally.

Turning finally to the last and least disturbed galaxy, NGC 1149, despite the bursty appearance of the SFH there’s no evidence for a major starburst in the cosmologically recent past. Whether an older starburst can be detected in this kind of modeling approach needs investigating.

One last set of graphs that may be useful. These show cumulative star formation histories — basically the cumulative sum of mass contributions starting from the oldest time bin. This is similar to a mass growth history which is a popular visualization. In my calculation of the latter the contributions are to the present day stellar mass, so an allowance for mass loss and remnant mass is made3these come from the source of the SSP models and are themselves models. Probably they are somewhat better than guesses. These things are basically black boxes to users.. The graphs are for the central regions only. Note the major virtue of these is that the contributions of major episodes of star formation can be estimated at a glance.

Cumulative star formation histories for central regions of 4 disturbed early type galaxies

To wrap up this part of the post 3 of these galaxies are compatible with the “modern merger hypothesis,” that is they have experienced centrally concentrated but spatially wide spread starbursts. The reason two of them don’t have post-starburst characteristics in the Hδ – D4000 plane is their starbursts are still underway. The current burst of star formation contributes about 5-10% of the mass in the central regions of these two. How much more is available is unknown (at least to me until I get around to finding out if there are HI mass estimates available).

Future plans: I’ve completed model runs on the 24 “post-starburst” galaxies in the MaNGA ancillary program dedicated to them. I may have something to say about them. I also may have something to say about one of the Zoogems targets that I had a small part in selecting.

Continue reading “Journal notes: Haines et al. (2015), “Testing the modern merger hypothesis…””

A little more on Schawinski’s blue early type galaxies

As I mentioned two posts ago there are 24 of these galaxies in the final MaNGA data release, a remarkable 11% of the full sample. I ran my SFH model code on all of these along with the prerequisite redshift offset routine1I actually completed these some time ago. I just haven’t had time to do much analysis or write about them. SDSS thumbnails of the sample are shown below. As expected none of these have significant spiral structure visible at SDSS resolution, but at least a few are noticeably disturbed.

thumbnails_blueetg
SDSS thumbnail images of Schawinski et al.’s blue early type galaxies in MaNGA final data release (SDSS DR17)

I’m just going to discuss a few topics in this post. I’ll save a more detailed discussion for when I’ve completed analysis of the ancillary post-starburst sample, which is underway now. First, here are velocity fields calculated for the stacked RSS data, with a signal to noise cutoff of 3, the same as I used for my analysis of rotation curves of disk galaxies. Note in the graphic below the ordering is different from the image thumbnails.

vfs_blueetg
Line of sight velocity fields of Schawinski et al.’s blue early type galaxies in the final MaNGA data release

By my count (based entirely on visual inspection) all but 2 of these exhibit large scale rotation, with perhaps 15 or 16 classifiable as regular rotators with the remainder containing multiple velocity components including a couple with (perhaps) kinematically distinct cores. The preponderance of rotating systems surprised me at first, but according to a review by Cappellari (2016) large scale rotation is predominant at least at lower stellar masses (Schawinski et al. characterized their sample as being “low to intermediate mass” among early type galaxies). The velocity fields indicate that many of these contain stellar disks, perhaps embedded in large bulges. That’s still consistent with classification as “early type galaxies.” Apparently the original Galaxy Zoo classification page used the term “elliptical” as the early type galaxy choice, but in the data release paper by Lintott et al. (2011) there’s a statement that the “elliptical” class should comprise ellipticals, S0’s, and perhaps Sa’s from Hubble’s classification scheme.

Depending on how my effort to do non-parametric line of sight velocity modeling goes I may return to examine the kinematics of this sample in more detail, in particular to look for evidence of gas and stellar kinematic decoupling.

Turning to the recent star formation history this sample runs the gamut from large scale starbursts to passively evolving as seen in the plot of (100 Myr averaged) star formation rate versus stellar mass density for all analyzed binned spectra (of which there were 1525 in the full sample). For reference the straight line is my estimate of the center of the local “spatially resolved star formation main sequence.” This is just a weighted least squares fit to the sample of 20 non-barred spirals with star forming BPT diagnostics that I discussed some time ago. My SFMS relation has the same slope as estimated by Bluck but is offset higher by about 0.7 dex, which probably just reflects the very different methods used to estimate star formation rates. The contour lines are the densest part of the relationship from the passively evolving Coma cluster sample that I also discussed in that post. The majority of the blue etg sample falls in the green valley, consistent with Schawinski et al.’s observation that only about 1/2 of the sample showed evidence for ongoing star formation.

sfr_mstar_blueetg
“Spatially resolved” star formation rate density versus stellar mass density for 24 blue early type galaxies in final MaNGA data release. Contour lines are corresponding values for 33 passively evolving Coma cluster galaxies.

Most of the points offset the most on the high side of the SFMS come from just two galaxies: MRK 888, which I’ve discussed in the last few posts, and SDSS J014143.18+134032.8 (this is apparently not in any “classical” catalog). The legacy survey cutout below clearly shows an extended tidal tail that’s a certain sign of a relatively recent merger.

SDSS J014143.18+134032.8, a disturbed, star-bursting blue early type galaxy

I just want to take a quick look at this one: below are maps of the star formation rate density and SSFR as well as scatterplots of the same against distance from the IFU center. As with MRK 888 ongoing star formation is widespread with a peak near the center, a classic case of a merger fueled starburst. In this galaxy star formation peaks in a ring somewhat outside the nucleus. The ring can be seen clearly in the SDSS cutout and must consist of HII regions.

8095-1902_sfr_ssfr
SDSS J014143.18+134032.8 (mangaid 1-41541; plateifu 8095-1902) Star formation rate density and specific star formation rate – maps and scatterplots against radius in kpc.

Schawinski et al. briefly discuss the possibility that their blue ETG’s could be progenitors of E+A (aka K+A) galaxies. This galaxy and MRK 888 are plausible candidates — if star formation shut off rapidly they would certainly exhibit strong Balmer absorption for a time after emission lines disappeared since they already do. Other members of this sample are already fading towards the red sequence, and if they ever qualified as “post-starburst” it must have been in the past.

I plan to look at star formation histories in more detail after I’ve completed model runs on the MaNGA post-starburst sample.

What fraction of Schawinski’s “Blue early type galaxies” are ellipticals?

The first iteration of Galaxy Zoo led to several collections of distinct objects, including a sample of 215 “blue early type galaxies” published in Schawinski et al. (2009)1which inexplicably and consistently says there were 204 objects while the catalog published in Vizier contains 215.I found this an interesting group of galaxies, partly because of a possible link to post-starburst (K+A) galaxies that was discussed in the original paper. The authors discuss at some length the likelihood that these are results of mergers in the cosmologically recent past, with at least one of the progenitors being gas rich. Many (at least 25% and possibly more than half) were found to be currently starforming and the rest likely to have only recently ceased forming stars as inferred from their blue colors.

The ongoing Zoogems program has 12 of Schawinski’s blue ETGs on its target list, of which 6 have been observed so far as of mid-January 2022. Somewhat surprisingly there are 24 in the final MaNGA release, over 11% of the sample!

Taking a look at the 6 with HST observations I would say none of these are typical ellipticals. Five show some degree of spiral structure although in 4 it’s embedded in a more diffuse body. One appears to me to be an S0 with both inner and outer rings — this is in agreement with the one published morphological classification I’ve found. All of the others appear more disky than ellipsoidal to me, although this is just my possibly flawed qualitative judgment. At least two are visibly disturbed. One (CGCG 315-014) is connected to a nearby galaxy with a long tidal tail as seen in the Legacy Survey thumbnail below. Markarian 888, which will be the subject of the rest of this post, has shells that extend well past the main body of the galaxy and prominent, centrally concentrated dust lanes.

CGCG 315-014 Legacy Survey Thumbnail

So far it’s the only Zoogems blue etg target with a MaNGA observation (two others on the target list are in MaNGA but of course there’s no guarantee they will ever be observed). As is often the case the IFU could have been larger — this was observed with a 37 fiber bundle giving 111 dithered spectra in the RSS file.

MRK 888 SDSS thumbnail with MaNGA IFU footprint

As always the first step in analyzing these data is to estimate redshift offsets for each spectrum, and from there we get a velocity field, which in this case shows a rapid rotator with a fairly symmetrical radial velocity pattern.

Mrk 888 (MaNGA plateifu 9894-3703) velocity field

Visual inspection suggests the line of sight velocity distribution is consistent with a rotating thin disk, so I fed the data to my Gaussian process based rotation modeling code, with results summarized below. In fact the model does an excellent job of accounting for the data, with residuals (not shown) from the model fit (top right) in a range of ±15 km/sec. One unusual feature of the velocity field is the rotation velocity turns over at somewhat less than one effective radius. Whether the rotation curve declines smoothly outside the IFU footprint or is kinematically disconnected from the outer parts of the galaxy is of course unknowable at this time.

Gaussian process rotation model results

I also ran my usual star formation history modeling code on the data binned to 97 spectra. First, here are some summary results. The stellar mass density declines roughly exponentially, which is consistent with a disky morphology:

Model estimate of stellar mass density vs. radius

Next are maps of the estimated Hα luminosity density and, on the right, the BPT classification from the [O III]/Hβ vs. [N II]/Hα diagnostic. The contours are elliptical with major axes closely aligned to the rotation axis (the posterior mean for the angle is the dashed line in the velocity field plot above). Again, the emission appears to arise in a disk.

The proper interpretation of the “composite” BPT classification is something I think I’ve written about in the past. It was originally suggested to indicate a mix of AGN and stellar ionization, but here it arises in a thin ring of weak but detectable emission just outside the star forming region. If it’s truly composite it’s likely to arise from a mix of weak star formation and ionization by hot evolved stars. In any case there’s no evidence for an AGN in the optical data.

(L) Hα luminosity density (R) BPT classification from [O III]/Hβ vs {N II]/Hα diagnostic

Next are maps of the modeled (100 Myr average) star formation rate density and specific star formation rate, and in the second row scatter plots of the same estimates against radius in kpc. The trends with radius are somewhat unusual, especially for SSFR which in a normal disk galaxy typically increases with radius even if the highest total star formation rates are centrally concentrated. Highly centrally concentrated star formation in the aftermath of mergers is predicted by some simulations.

(TL) star formation rate density; (TR) specific star formation rate; (bottom row) scatter plots vs. radius

A couple more graphs will round out my discussion of summary model estimates. As I’ve shown several times before there’s a pretty tight linear relationship between modeled SFR density and estimated Hα luminosity density. In this plot Hα is corrected for modeled stellar attenuation, which is expected always to underestimate the attenuation in emission line emitting regions. That, and the fact that Hα emission and the model star formation rate estimates probe order of magnitude different time scales probably account for the systematic offset from the standard calibration given by the straight line.

Model star formation rate density vs. Hα luminosity density corrected for stellar attenuation. Straight line is calibration from Calzetti (2012).

And, once again I show a map of the modeled optical depth of stellar attenuation. The region of highest optical depth nicely tracks the visible dust (the HST image at the top is rotated about 90º from the SDSS image). Outside the dusty region there appears to be a shallow gradient, which might indicate that the nearer side is to the northeast.

Map of modeled optical depth of stellar attenuation

Finally here are plots of the model star formation history for all spectra ordered by distance from the IFU center. In the inner 1.5 kpc or so there’s some recent burstiness with possibly a very recent acceleration of star formation. For reasons I’ve discussed recently I don’t take either the timing or magnitude of bursts of star formation too seriously, but the behavior of the models is consistent with a recent revival of star formation due presumably to a merger, for which there are multiple lines of evidence.

model star fomation histories for all spectra

With 24 of these galaxies and another 31 from the compilation of Melnick and dePropris and the post-starburst ancillary program in the final release of MaNGA these samples satisfy my criteria of being manageably sized for my computing resources while large enough to say something about the groups. So, when time permits I plan to take a look. I already have the data in hand.

KUG 0859+406 – unravelling the differences between 2018 and 2021 model runs

I did have my old data and model runs of course, in fact they were spread over several directories on two machines. I’m going to refer to it by this catalog designation, KUG standing for the “Kiso survey of Ultraviolet-excess Galaxies.” It’s also a low power radio source with catalog entries in both FIRST and NVSS, and of course it’s in MaNGA with plateifu ID 8440-6104 (mangaid 01-216976).

In my 2018 model runs, which were interesting enough to write 3 posts about, I found this galaxy had undergone an extraordinarily large burst of star formation that began ~1 Gyr ago with locally as much as 60% of the present day stellar mass born in the burst and something like 40% of the mass over the footprint of the IFU. In this years model runs the peak burst fraction was a considerably more modest ~25% and globally barely amounted to a slight enhancement of star formation. The starburst was also much more localized than in the earlier runs:

burst_strengths_compared
Fractional stellar mass in stars between 0.1 and 1.75 Gyr old in 2018 and 2021 model runs

So what happened? First, here is a comparison of modeled star formation histories for the innermost fiber, which got the largest injection of mass in the starburst.

central_sfh_compared
Model star formation histories for central fiber of MaNGA plateifu 8440-6104, 2018 and 2021 model runs

The obvious remark is the double peaked starburst noted back in 2018 (and discussed at some length) has been replaced with a single narrow peak with a slow ramp up and fast decay. The peak SFR is a little larger than before but the total mass in the burst is lower.

I’ve made several changes in model formulation since 2018, of which the most important in the current context is adopting the more flexible “modified Calzetti” attenuation relation that adds an additional slope parameter to the prescription. In the current year model runs a steeper than Calzetti relation is favored throughout the IFU footprint, particularly in the central region where the starburst was strongest:

map_delta
Map of modified Cal;zetti slope parameter δ — MaNGA plateifu 8440-6104

A smaller optical depth of attenuation is also favored throughout:

tauv_compared
Modeled optical depth of attenuation – 2021 runs vs. 2018 MaNGA plateifu 8440-6104

This has a couple predictable consequences. Steeper attenuation will favor an intrinsically bluer, hence younger population while a lower optical depth requires less light, and hence mass in the stellar population. I can test this directly by returning to a model with Calzetti attenuation, and here is the result for the central fiber (this model run is labeled 2021 (c) in the legend below):

mgh_compared
Mass growth histories – 2021 run 2021 run with Calzetti attenuation 2018 run Central fiber of MaNGA plateifu 8440-6104

So, an eyeball analysis suggests about 3/4 of the difference between the 2018 and 2021 runs is due to the modification to the attenuation relation. The other changes I’ve made to the models are to change the stellar contribution parameters from a non-negative vector to a simplex, and at the same time changing the way I rescale the data. In early runs the SSP model fluxes were scaled to make the maximum stellar contribution ≈ 1, while the current models scale both the galaxy and SSP fluxes to ≈ 1 in the neighborhood of V, making the individual stellar contributions approximately the fraction of light contributed. An additional scale factor parameter in the model is used to adjust the overall fit. Assuming I did this right this should have no effect on a deterministic maximum likelihood solution, but with MCMC who knows?

Although the fit to the data looks about the same between the model with and without the attenuation modification the summed log-likelihood is consistently about 1% higher for the modified Calzetti model with no overlap at all in the distribution of likelihood. This suggests the case for a steeper than Calzetti attenuation is a fairly robust result.

ppfits
“Posterior predictive” fits to galaxy flux data – modified Calzetti attenuation vs. Calzetti – central fiber of MaNGA plateifu 8440-6104

The galaxy flux data also changed a little bit. The early runs were on the DR14 release (version 2_1_2 of the MaNGA DRP) while the recent ones used the DR15 release (ver 2_4_3). Most of the calibration differences resemble random noise, but there is some curvature that systematically affects both the red and blue ends of the spectrum and could cause some change in the temperature distribution of the models:

measured_flux_compared
Difference in measured flux from DR14 to DR15 – central fiber of MaNGA plateifu 8440-6104

While the detailed star formation histories changed, quantities that aren’t too model dependent didn’t very much. One example is shown below. Also, the kinematic maps agree with the earlier ones in detail.

halpha_luminosity_compared
Hα luminosity density – 2021 runs vs. 2018 MaNGA plateifu 8440-6104
velocity_field_8440_6104
Velocity field of MaNGA plateifu 8440-6104 from 2021 model runs. Map interpolated from RSS file spectra.

One input that hasn’t changed are the emiles SSP model spectra, although there have been some procedural changes in how I handle the modeling. Early on I often used a much smaller subset of SSP models with just 27 time bins and 3 metallicities for preliminary modeling, including my first models on the same binning of these data. I also routinely ran 250 warmup iterations with just 250 more per chain. My current standard practice is always to use the largest emiles subset with 216 SSP models in 54 time bins and 4 metallicities, and I generally run 750 post-warmup iterations per chain but still with only 250 warmup iterations. This is generally enough and if adaptation fails it is usually fairly obvious. The small sample size of the earlier runs mostly effects the precision of inferences rather than means.

To conclude for now, my speculation about whether it might be possible to say something about the timing of critical events in a merger from the model star formation history was too optimistic. On a positive note though both sets of model runs retrodict that coalescence occurred at a lookback time around 500Myr ago, which is consistent with the fact that tidal tails and other merger signatures are clearly visible even in SDSS imaging. Both sets of model runs also have that odd uptick in star formation at 30Myr in the central fiber. And while the difference in burst mass contributions is a little disconcerting the current runs are more consistent with the likely gas content of ordinary spiral galaxies.

This example illustrates another well known “degeneracy” among attenuation (and adopted attenuation relation), mass, and stellar age. Whether I’ve broken the degeneracy by adopting the more flexible attenuation prescription described some time ago remains to be validated.

Markarian 848 – Closing topics

I’m going to close out my analysis of Mrk 848 for now with three topics. First, dust. Like most SED fitting codes mine produces an estimate of the internal attenuation, which I parameterize with τV, the optical depth at V assuming a conventional Calzetti attenuation curve. Before getting into a discussion for context here is a map of the posterior mean estimate for the higher S/N target binning of the data. For reference isophotes of the synthesized r band surface brightness taken from the MaNGA data cube are superimposed:

mrk848_tauv_map
Map of posterior mean of τV from single component dust model fits with Calzetti attenuation

This compares reasonably well with my visual impression of the dust distribution. Both nuclei have very large dust optical depths with a gradual decline outside, while the northern tidal tail has relatively little attenuation.

The paper by Yuan et al. that I looked at last time devoted most of its space to different ways of modeling dust attenuation, ultimately concluding that a two component dust model of the sort advocated by Charlot and Fall (2000) was needed to bring results of full spectral fitting using ppxf on the same MaNGA data as I’ve examined into reasonable agreement with broad band UV-IR SED fits.

There’s certainly some evidence in support of this. Here is a plot I’ve shown for some other systems of the estimated optical depth of the Balmer emission line emitting regions based on the observed vs. theoretical Balmer decrement (I’ve assumed an intrinsic Hα/Hβ ratio of 2.86 and a Calzetti attenuation relation) plotted against the optical depth estimated from the SFH models, which roughly estimates the amount of reddening needed to fit the SSP model spectra to the observed continuum. In some respects this is a more favorable system than some I’ve looked at because Hβ emission is at measurable levels throughout. On the other hand there is clear evidence that multiple ionization mechanisms are at work, so the assumption of a single canonical value of Hα/Hβ is likely too simple. This might be a partial cause of the scatter in the modeled relationship, but it’s encouraging that there is a strong positive correlation (for those who care, the correlation coefficient between the mean values is 0.8).

The solid line in the graph below is 1:1. The semi-transparent cloud of lines are the sampled relationships from a Bayesian errors in variables regression model. The mean (and marginal 1σ uncertainty) is \(\tau_{V, bd} = (0.94\pm 0.11) + (1.21 \pm 0.12) \tau_V\). So the estimated relationship is just a little steeper than 1:1 but with an offset of about 1, which is a little different from the Charlot & Fall model and from what Yuan et al. found, where the youngest stellar component has about 2-3 times the dust attenuation as the older stellar population. I’ve seen a similar not so steep relationship in every system I’ve looked at and don’t know why it differs from what is typically assumed. I may look into it some day.

τV estimated from Balmer decrement vs. τV from model fits. Straight line is 1:1 relation. Cloud of lines are from Bayesian errors in variables regression model.

I did have time to run some 2 dust component SFH models. This is a very simple extension of the single component models: a single optical depth is applied to all SSP spectra. A second component with the optical depth fixed at 1 larger than the bulk value is applied only to the youngest model spectra, which recall were taken from unevolved SSPs from the updated BC03 library. I’m just going to show the most basic results from the models for now in the form of maps of the SFR density and specific star formation rate. Compared to the same maps displayed at the end of the last post there is very little difference in spatial variation of these quantities. The main effect of adding more reddened young populations to the model is to replace some of the older light — this is the well known dust-age degeneracy. The average effect was to increase the stellar mass density (by ≈ 0.05 dex overall) while slightly decreasing the 100Myr average SFR (by ≈ 0.04 dex), leading to an average decrease in specific star formation rate of ≈ 0.09 dex. While there are some spatial variations in all of these quantities no qualitative conclusion changes very much.

mrk848_sigma_sfr_sfr_2dust_maps
Results from fits with 2 component dust models. (L) SFR density. (R) Specific SFR

Contrary to Yuan+ I don’t find a clear need for a 2 component dust model. Without trying to replicate their results I can’t say why exactly we disagree, but I think they erred in aggregating the MaNGA data to the lowest spatial resolution of the broad band photometric data they used, which was 5″ radius. There are clear variations in physical conditions on much smaller scales than this.

Second topic: the most widely accepted SFR indicator in visual wavelengths is Hα luminosity. Here is another plot I’ve displayed previously: a simple scatterplot of Hα luminosity density against the 100Myr averaged star formation rate density from the SFH models. Luminosity density is corrected for attenuation estimated from the Balmer decrement and for comparison the light gray points are the uncorrected values. Points are color coded by BPT class determined in the usual manner. The straight line is the Hα – SFR calibration of Moustakas et al. (2006), which in turn is taken from earlier work by Kennicutt.

Model SFR density vs. Hα luminosity density corrected for extinction estimated from Balmer decrement. Light colored points are uncorrected for extinction. Straight line is Hα-SFR calibration from Moustakas et al. (2006)

Keeping in mind that Hα emission tracks star formation on timescales of ∼10 Myr1to the extent that ionization is caused by hot young stars. There are evidently multiple ionizing sources in this system, but disentangling their effects seems hard. Note there’s no clear stratification by BPT class in this plot. this graph strongly supports the scenario I outlined in the last post. At the highest Hα luminosities the SFR-Hα trend nicely straddles the Kennicutt-Moustakas calibration, consistent with the finding that the central regions of the two galaxies have had ∼constant or slowly rising star formation rates in the recent past. At lower Hα luminosities the 100Myr average trends consistently above the calibration line, implying a recent fading of star formation.

The maps below add some detail, and here the perceptual uniformity of the viridis color palette really helps. If star formation exactly tracked Hα luminosity these two maps would look the same. Instead the northern tidal tail in particular and the small piece of the southern one within the IFU footprint are underluminous in Hα, again implying a recent fading of star formation in the peripheries.

(L) Hα luminosity density, corrected for extinction estimated by Balmer decrement. (R) SFR density (100 Myr average).

Final topic: the fit to the data, and in particular the emission lines. As I’ve mentioned previously I fit the stellar contribution and emission lines simultaneously, generally assuming separate single component gaussian velocity dispersions and a common system velocity offset. This works well for most galaxies, but for active galaxies or systems like this one with complex velocity profiles maybe not so much. In particular the northern nuclear region is known to have high velocity outflows in both ionized and neutral gas due presumably to supernova driven winds. I’m just going to look at the central fiber spectrum for now. I haven’t examined the fits in detail, but in general they get better outside the immediate region of the center. First, here is the fit to the data using my standard model. In the top panel the gray line, which mostly can’t be seen, is the observed spectrum. Blue are quantiles of the posterior mean fit — this is actually a ribbon, although its width is too thin to be discernable. The bottom panel are the residuals in standard deviations. Yes, they run as large as ±50σ, with conspicuous problems around all emission lines. There are also a number of usually weak emission lines that I don’t track that are present in this spectrum.

mrk848_fit_central_spec
Fit to central fiber spectrum; model with single gaussian velocity distributions.

I have a solution for cases like this which I call partially parametric. I assume the usual Gauss-Hermite form for the emission lines (as in, for example, ppxf) while the stellar velocity distribution is modeled with a convolution kernel2I think I’ve discussed this previously but I’m too lazy to check right now. If I haven’t I’ll post about it someday. Unfortunately the Stan implementation of this model takes at least an order of magnitude longer to execute than my standard one, which makes its wholesale use prohibitively expensive. It does materially improve the fit to this spectrum although there are still problems with the stronger emission lines. Let’s zoom in on a few crucial regions of the spectrum:

Zoomed in fit to central fiber spectrum using “partially parametric velocity distribution” model. Grey: observed flux. Blue: model.

The two things that are evident here are the clear sign of outflow in the forbidden emission lines, particularly [O III] and [N II], while the Balmer lines are relatively more symmetrical as are the [S II] doublet at 6717, 6730Å. The existence of rather significant residuals is likely because emission is coming from at least two physically distinct regions while the fit to the data is mostly driven by Hα, which as usual is the strongest emission line. The fit captures the emission line cores in the high order Balmer lines rather well and also the absorption lines on the blue side of the 4000Å break except for the region around the [Ne III] line at 3869Å.

I’m mostly interested in star formation histories, and it’s important to see what differences are present. Here is a comparison of three models: my standard one, the two dust component model, and the partially parametric velocity dispersion model:

mrk848_centralsfr3ways
Detailed star formation history models for the northern nucleus using 3 different models.

In fact the differences are small and not clearly outside the range of normal MCMC variability. The two dust model slightly increases the contribution of the youngest stellar component at the expense of slightly older contributors. All three have the presumably artifactual uptick in SFR at 4Gyr and very similar estimated histories for ages >1 Gyr.

I still envision a number of future model refinements. The current version of the official data analysis pipeline tracks several more emission lines than I do at present and has updated wavelengths that may be more accurate than the ones from the SDSS spectro pipeline. It might be useful to allow at least two distinct emission line velocity distributions, with for example one for recombination lines and one for forbidden. Unfortunately the computational expense of this sort of generalization at present is prohibitive.

I’m not impressed with the two dust model that I tried, but there may still be useful refinements to the attenuation model to be made. A more flexible form of the Calzetti relation might be useful for example3there is recent relevant literature on this topic that I’m again too lazy to look up.

My initial impression of this system was that it was a clear false positive that was selected mostly because of a spurious BPT classification. On further reflection with MaNGA data available it’s not so clear. A slight surprise is the strong Balmer absorption virtually throughout the system with evidence for a recent shut down of star formation in the tidal tails. A popular scenario for the formation of K+A galaxies through major mergers is that they experience a centrally concentrated starburst after coalescence which, once the dust clears and assuming that some feedback mechanism shuts off star formation leads to a period of up to a Gyr or so with a classic K+A signature4I previously cited Bekki et al. 2005, who examine this scenario in considerable detail.Capturing a merger in the instant before final coalescence provides important clues about this process.

To the best of my knowledge there have been no attempts at dynamical modeling of this particular system. There is now reasonably good kinematic information for the region covered by the MaNGA IFU, and there is good photometric data from both HST and several imaging surveys. Together these make detailed dynamical modeling technically feasible. It would be interesting if star formation histories could further constrain such models. Being aware of the multiple “degeneracies” between stellar age and other physical properties I’m not highly confident, but it seems provocative that we can apparently identify distinct stages in the evolutionary history of this system.

Markarian 848 – detailed star formation histories

In this post I’m going to take a quick look at detailed, spatially resolved star formation histories for this galaxy pair and briefly compare to some of the recent literature. Before I start though, here is a reprocessed version of the blog’s cover picture with a different crop and some Photoshop curve and level adjustments to brighten the tidal tails a bit. Also I cleaned up some more of the cosmic ray hits.

mrk_848_hst_crop_square
Markarian 848 – full resolution crop with level curve adjustment HST ACS/WFC F435W/F814W/F160 false color image

The SFH models I’m going to discuss were based on MaNGA data binned to a more conservative signal to noise target than in the last post. I set a target S/N of 25 — Cappellari’s Voronoi binning algorithm has a hard coded acceptance threshold of 80% of the target S/N, which results in all but a few bins having an actual S/N of at least 20. This produced a binned data set with 63 bins. The map below plots the modeled (posterior mean log) stellar mass density, with clear local peaks in the bins covering the positions of the two nuclei. The bins are numbered in order of increasing distance of the weighted centroids from the IFU center, which corresponds with the position of the northern nucleus. The center of bin 1 is slightly offset from the northern nucleus by about 2/3″. For reference the angular scale at the system redshift (z ≈ 0.040) is 0.8 kpc/” and the area covered by each fiber is 2 kpc2.

Although there’s nothing in the Voronoi binning algorithm that guarantees this it did a nice job of partitioning the data into physically distinct regions, with the two nuclei, the area around the northern nucleus, and the bridge between them all binned to single fiber spectra. The tidal tails are sliced into several pieces while the very low surface brightness regions to the NE and SW get a few bins.

mrk848_stmass_map
Stellar mass density and bin numbers ordered in increasing distance from northern nucleus

I modeled star formation histories with my largest subset of the EMILES library consisting of SSP model spectra for 54 time bins and 4 metallicity bins. As in previous exercises I fit emission lines and the stellar contribution simultaneously, which is riskier than usual in this system because some regions, especially near the northern nucleus, have complex velocity structures producing emission line profiles that are far from single component gaussians. I’ll look at this in more detail in a future post, but for now let’s just take these SFH models at face value and see where they lead. As I’ve done in several previous posts what’s plotted below are star formation rates in solar masses/year over cosmic time, counting backwards from “now” (now being when the light left the galaxies about 500Myr ago). This time I’ve used log scaling for both axes and for the first time in these posts I’ve displayed results for every bin ordered by increasing distance from the IFU center. The first two rows cover the central ~few kpc region surrounding the northern nucleus. The southern nucleus covered by bin 44 is in row 8, second from left. Both the time and SFR scales are the same for all plots. The black lines are median posterior marginal estimates and the ribbons are 95% posterior confidence intervals.

mrk848_sfh_bin20
Star formation histories for the Voronoi binned regions shown in the previous map. Numbers at the top correspond to the bin numbers in the map.

Browsing through these, all regions show ∼constant or gradually decreasing star formation up to about 1 Gyr ago1You may recall I’ve previously noted there is always an uptick in star formation rate at 4 Gyr in my EMILES based models, and that’s seen in all of these as well. This must be spurious, but I still don’t know the cause.. This of course is completely normal for spiral galaxies.

Most regions covered by the IFU began to show accelerated and in some areas bursts of star formation at ∼1 Gyr. In more than half of the bins the maximum star formation rate occurred around 40-60 Myr ago, with a decline or in some areas cessation of star formation more recently. In the central few kpc around the northern nucleus on the other hand star formation accelerated rapidly beginning ~100 Myr ago and continues at high levels. The southern nucleus has a rather different estimated recent star formation history, with no (visible) starburst and instead gradually increasing SFR to a recent maximum. Ongoing star formation at measurable levels is much more localized in the southern central regions and weaker by factors of several than the northern companion.

Here’s a map that I thought would be too noisy to be informative, but turns out to be rather interesting. This shows the spatial distribution of the lookback time to the epoch of maximum star formation rate estimated by the marginal posterior mean. The units displayed are log(age) in Myr; recall that I added unevolved SSP models from the 2013 update of the BC03 models to the BaSTI based SSP library, assigning them an age of 10Myr, so a value of 1 here basically means ≈ now.

mrk848_lbt_to_maxsfr
Look back time to epoch of maximum star formation rate as estimated by marginal posterior mean. Units are log(age in Myr).

To summarize, there have been three phases in the star formation history of this system: a long period of normal disk galaxy evolution; next beginning about 1 Gyr ago widespread acceleration of the SFR with localized bursts; and now, within the past 50-100 Myr a rapid increase in star formation that’s centrally concentrated in the two nuclei (but mostly the northern one) while the peripheral regions have had suppressed activity.

I previously reviewed some of the recent literature on numerical studies of galaxy mergers. The high resolution simulations of Hopkins et al. (2013) and the movies available online at http://www.tapir.caltech.edu/~phopkins/Site/Movies_sbw_mgr.html seem especially relevant, particularly their Milky Way and Sbc analogs. They predict a general but clumpy enhancement of star formation at around the time of first pericenter passage; a stronger, centrally concentrated starburst at the time of second pericenter passage, with a shorter period of separation before coalescence. Surprisingly perhaps, the merger timescales for both their MW and Sbc analogs are very similar to my SFH models, with ∼ 1 Gyr between first and second perigalactic passage and another few 10’s of Myr to final coalescence.

I’m going to wrap up with maps of the (posterior mean estimates) of star formation rate density (as \(\log10(\mathsf{M_\odot/yr/kpc^2})\)) and specific star formation rate, which has units log10(yr-1). These are 100Myr averages and recall are based solely on the SSP template contributions.

mrk848_sfr_ssfr
(L) Star formation rate density (R) Specific star formation rate

Several recent papers have made quantitative estimates of star formation rates in this system, which I’ll briefly review. Detailed comparisons are somewhat difficult because both the timescales probed by different SFR calibrators differ and the spatial extent considered in the studies varies, so I’ll just summarize reported results and compare as best as I can.

Yuan et al. (2018) modeled broad band UV-IR SED’s taking data from GALEX, SDSS, and Spitzer; and also did full spectrum fitting using ppxf on the same MaNGA data I’ve been examining. They divided the data into 5″ (≈ 4 kpc) radius regions covering the tidal tails and each nucleus. The broadband data were fit with parametric SFH models with a pair of exponential bursts (one fixed at 13Gyr, the other allowed to vary). From the parametric models they estimated the SFR in the northern and southern nuclei as 76 and 11 \(\mathsf{M_\odot/yr}\) (the exact interpretation of this value is unclear to me). For comparison I get values of 33 and 13 \(\mathsf{M_\odot/yr}\) for regions within 4 kpc of the two nuclei by calculating the average (posterior mean) star forming density and multiplying by π*42. They also calculated 100Myr average SFRs from full spectrum fitting with several different dust models and also from the Hα luminosity, deriving estimates that varied by an order of magnitude or more. Qualitatively we reach similar conclusions that the tails had earlier starbursts and are now forming stars at lower rates than their peak, and also that the northern nucleus has recent star formation a factor of several higher than the southern.

Cluver et al. (2017) were primarily trying to derive SFR calibrations for WISE W3/W4 filters using samples of normal star forming and LIRG/ULIRGs (including this system) with Spitzer and Herschel observations. Oddly, although they tabulate IR luminosities for the entire sample they don’t tabulate SFR estimates. But plugging into their equations 5 and 7 I get star formation rate estimates of just about 100 \(\mathsf{M_\odot/yr}\). These are global values for the entire system. For comparison I get a summed SFR for my models of ≈ 45 \(\mathsf{M_\odot/yr}\) (after making a 0.2 dex adjustment for fiber overlap).

Tsai and Hwang (2015) also used IR data from Spitzer and conveniently present results for quantities I track using the same units. Their estimates of the (log) stellar mass density, SFR density, and specific star formation rate in the central kpc of (presumably) the northern nucleus were 9.85±0.09, 0.43±0.00, and -9.39±0.09 respectively. For comparison in the fiber that covers the northern nucleus my models return 9.41±0.02, 0.57±0.04, and -8.84±0.04. For the 6 fibers closest to the nucleus the average stellar mass density drops to 9.17 and SFR density to about 0.39. So, our SFR estimates are surprisingly close while my mass density estimate is as much as 0.7 dex lower.

Finally, Vega et al. (2008) performed starburst model fits to broadband IR-radio data. They estimated a burst age of about 60Myr for this system, with average SFR over the burst of 225 \(\mathsf{M_\odot/yr}\) and current (last 10 Myr) SFR 87 \(\mathsf{M_\odot/yr}\) in a starforming region of radius 0.27kpc. Their model has optical depth of 33 at 1 μm, which would make their putative starburst completely invisible at optical wavelengths. Their calculations assumed a Salpeter IMF, which would add ≈0.2 dex to stellar masses and star formation rates compared to the Kroupa IMF used in my models.

Overall I find it encouraging that my model SFR estimates are no worse than factors of several lower than what are obtained from IR data — if the Vega+ estimate of the dust optical depth is correct most of the star formation is well hidden. Next time I plan to look at Hα emission and tie up other loose ends. If time permits before I have to drop this for a while I may look at two component dust models.

Markarian 848

This galaxy1also known as VV705, IZw 107, IRAS F15163+4255, among others has been my feature image since I started this blog. Why is that besides that it’s kind of cool looking? As I’ve mentioned before I took a shot a few years ago at selecting a sample of “transitional” galaxies from the SDSS spectroscopic sample, that is ones that may be in the process of rapidly shutting down star formation2See for example Alatalo et al. (2017) for recent usage of this terminology.. I based the selection on a combination of strong Hδ absorption and either weak emission lines or line ratios other than pure starforming, using measurements from the MPA-JHU pipeline. This galaxy pair made the sample based on the spectrum centered on the northern nucleus (there is also a spectrum of the southern nucleus from the BOSS spectrograph, but no MPA pipeline measurements). Well now, these galaxies are certainly transitioning to something, but they’re probably not shutting down star formation just yet. Simulations of gas rich mergers generally predict a starburst that peaks around the time of final coalescence. There is also significant current star formation, as high as 100-200 \(M_\odot/yr\) per various literature estimates, although it is mostly hidden. So on the face of it at least this appears to be a false positive. This was also an early MaNGA target, and one of a small number where both nuclei of an ongoing merger are covered by a single IFU:

Markarian 848. SDSS image thumbnail with MaNGA IFU footprint.

Today I’m going to look at a few results of my analysis of the IFU data that aren’t too model dependent to get some insight into why this system was selected. As usual I’m looking at the stacked RSS data rather than the data cubes, and for this exercise I Voronoi binned the spectra to a very low target SNR of 6.25. This leaves most of the fibers unbinned in the vicinity of the two nuclei. First, here is a map of BPT classification based on the [O III]/Hβ vs. [N II]/Hα diagnostic as well as scatterplots of several BPT diagnostics. Unfortunately the software I use for visualizing binned maps extends the bins out to the edge of the bounding rectangle rather than the hexagonally shaped edge of the data, which is outlined in red. Also note that different color coding is used for the scatter plots than the map. The contour lines in the map are arbitrarily spaced isophotes of the synthesized R band image supplied with the data cube.

Map of BPT class determined from [O III]/Hβ vs. [N II}/Hα diagnostic diagram, and BPT diagnostics for [N II], [S II] and [O I]. Curves are boundary lines form Kewley et al. (2006).

What’s immediately obvious is that most of the area covered by the IFU including both nuclei fall in the so-called “composite” region of the [N II]/Hα BPT diagram. This gets me back to something I’ve complained about previously. There was never a clear physical justification for the composite designation (which recall was first proposed in Kauffmann et al. 2003), and the upper demarcation line between “pure” AGN and composite systems as shown in the graph at upper right was especially questionable. It’s now known if it wasn’t at the turn of the century (which I doubt is the case) that a number of ionization mechanisms can produce line ratios that fall generally in the composite/LINER regions of BPT diagnostics. Shocks in particular are important in ongoing mergers. High velocity outflow of both ionized and neutral gas have been observed in the northern nucleus by Rupke and Veillux (2013), which they attributed to supernova driven winds.

The evidence for AGN in either nucleus is somewhat ambiguous. Fu et al. (2018) called this system a binary AGN, but that was based on the “composite” BPT line ratios from the same MaNGA data as we are examining (their map, by the way, is nearly identical to mine; see also Yuan et al. 2018). By contrast Vega et al. 2008 were able to fit the entire SED from NIR to radio frequencies with a pure starburst model and no AGN contribution at all, while more recently Dietrich et al. 2018 estimated the AGN fraction to be around 0.25 from NIR to FIR SED fitting. A similar conclusion that both nuclei contain both a starburst and AGN component was reached by Vardoulaki et al. 2015 based on radio and NIR data. One thing I haven’t seen commented on in the MaNGA data that possibly supports the idea that the southern nucleus harbors an AGN is that the regions with unambiguous AGN-like optical emission line ratios are fairly symmetrically located on either side of the southern nucleus and separated from it by ∼1-2 kpc. This could perhaps indicate the AGN is obscured from our view but visible from other angles.

There are also several areas with starforming emission line ratios just to the north and east of the northern nucleus and scattered along the northern tidal tail (the southern tail is largely outside the IFU footprint). In the cutout below taken from a false color composite of the F435W+F814W ACS images several bright star clusters can be seen just outside the more heavily obscured nuclear region, and these are likely sources of ionizing photons.

mrk_848_hst_crop
Markarian 848 nuclei. Cutout from HST ACS F814W+F435W images.

Finally turning to the other component of the selection criteria, here is a map of the (pseudo) Lick HδA index and a plot of the widely used HδA vs Dn(4000) diagnostic. It’s a little hard to see a clear pattern in the map because this is a rather noisy index, but strong Balmer absorption is seen pretty much throughout, with the highest values outside the two nuclei and especially along the northern tidal tail.

Location in the Hδ – D4000 plane doesn’t uniquely constrain the star formation history, but the contour plot taken from a largish sample of SDSS spectra from the NGP is clearly bimodal, with mostly starforming galaxies at upper left and passively evolving ones at lower right, with a long “green valley” in between. Simple models of post-starburst galaxies will loop upwards and to the right in this plane as the starburst ages before fading towards the green valley. This is exactly where most of points in this diagram lie, which certainly suggests an interesting recent star formation history.

mrk848_hda_d4000
(L) Map of HδA index. (R) HδA vs Dn(4000) index. Contours are for a sample of SDSS spectra from the north galactic pole.

I’m going to end with a bit of speculation. In simulations of gas rich major mergers the progenitors generally complete 2 or 3 orbits before final coalescence, with some enhancement of star formation during the perigalactic passages and perhaps some ebbing in between. This process plays out over hundreds of Myr to some Gyr. What I think we are seeing now is the 2nd or third encounter of this pair, with the previous encounters having left an imprint in the star formation history.

I’ve done SFH modeling for this binning of the data, and also for data binned to higher SNR and modeled with higher time resolution models. Next post I’ll look at these in more detail.