ASTEROID LIGHT CURVES FROM THE PALOMAR TRANSIENT FACTORY SURVEY: ROTATION PERIODS AND PHASE FUNCTIONS FROM SPARSE PHOTOMETRY

Building upon the work of Kaasalainen et al. (2001), Hanuš & Ďurech (2012) discuss the inversion of asteroid light-curve data taken over several oppositions to obtain a 3D shape solution. The form of δ (see Equation (1)) in this case consists of a large number of free parameters (several tens to hundreds). Results from inversion agree well with those from stellar occultations, adaptive optics imaging, and in situ spacecraft imagery (Hanuš et al. 2013). Knowledge of the detailed irregular shapes of asteroids improves our ability to constrain models of their internal structure, as well as the magnitude and timescale of spin and orbital evolution due to solar-radiation and thermal emission, including the Yarkovsky and YORP effects (see Bottke et al. 2013 and references therein).

A simpler model for δ—suitable for fitting to data sparser than those required for most inversion methods—is a Jacobi ellipsoid (Chandrasekhar 1969) in its principal-axis spin state. The light curve of such an ellipsoid is a double-peaked sinusoid, given by a simple expression depending solely (assuming constant surface albedo) on the axes ratio and angle between the line of sight and spin axis. The fitted amplitude thus yields a lower-bound elongation estimate for the asteroid.

The predicted distribution of the rotation frequencies of a collisionally equilibrated system of particles has long been claimed to be a Maxwellian function (Salo 1987), which—as reviewed by Pravec et al. (2002)—very well approximates the observed distribution of several hundred of the brightest (∼40 km or larger) asteroids, but breaks down for smaller objects, among which an excess of slow and fast rotators appear to exist. Steinberg & Sari (2015) more recently argue that collision instead leads to a Lévy distribution, and that a significant primordial spin component remains in the present observed population. Some studies that have examined the spin distribution of small objects are Pravec et al. (2008), Polishook & Borsch (2009), the Thousand Asteroid Light Curve Survey (Masiero et al. 2009), and two brief observing runs conducted within the Palomar Transient Factory⁷ (PTF) survey (Polishook et al. 2012; Chang et al. 2014a).

Warner et al. (2009) describe the Light Curve Database (LCDB), which compiles several thousand densely sampled light curves of asteroids targeted by dedicated observing teams. Light curves in the LCDB have the following features:

• 1.  
  LCDB light curves' dense sampling generally permits fitting of Fourier series with many harmonic terms;
• 2.  
  LCDB light curves are often sampled over the shortest time window necessary to measure the period and therefore generally do not require large or uncertain corrections due to phase-angle effects;
• 3.  
  LCDB light curves' fitted periods are assigned integer quality codes by a human reviewer (from 1 = poor to 3 = confident).

All three of the above features are either impractical or infeasible when the set of light curves is very large and the data sparsely sampled, as is the case for PTF. In this work we adopt the following modified approaches when fitting light curves:

• 1.  
  We truncate the rotation curve's Fourier-series fit after the 2nd harmonic, a simplification broadly justified by Harris et al. (2014) and the assumption of an ellipsoidal shape (see Section 3.1.2).
• 2.  
  We simultaneously fit a phase-function model with the rotational part.
• 3.  
  We use a machine-learned classifier to objectively aid in estimating the validity of each fitted period. The classifier is trained using all fitted light curves that have previously (and confidently) measured LCDB periods and takes into account the accuracy with which the true period was retrieved along with 20 light-curve metrics (fitted period, amplitude, ratio of peaks, per degree of freedom of fit, number of data points, and more).

Use of a machine classifier in asteroid light-curve period quality assessment is entirely novel and inspired in part by work done by PTF collaborators in extragalactic transient science (Bloom et al. 2012) and variable star science (Masci et al. 2014; Miller et al. 2014), as well as work by A. Waszczak et al. (2015, in preparation) on detection techniques for streaking near-Earth objects. Among the advantages of using a machine-classified quality score is that, via cross-validation with the known-period sample, one estimates the completeness and contamination, i.e., the true-positive rate (TPR) and false-positive rate (FPR) with respect to identifying an accurately fit period, as a function of, e.g., the period, amplitude, etc. The resulting TPR and FPR may then be used to de-bias the classifier-filtered period distribution.

The analytic phase function of an ideal Lambertian-scattering sphere fits well to featureless, atmospheric planets like Venus, but quite poorly to airless bodies (see Figure 3.9 of Seager 2010 for a comparison). In later sections we describe several ϕ models that have been derived for (or empirically fit to) asteroids. Qualitatively, asteroids show an approximately linearly decreasing ϕ out to , modified by a surge (increase in slope) at low phase angles (), known as the opposition effect (see Figure 1).

Zoom In Zoom Out Reset image size

Early work (e.g., Bowell et al. 1989 and references therein) on a small sample of well-observed asteroids suggested that different asteroid spectral types display distinct behavior in ϕ. Figure 1 compares example phase-curve data for D, C, S, and E types,⁸ incorporating photometry from various sources. We emphasize the fact that all of the data points in Figure 1 have been corrected for rotational modulation (the δ in Equation (1)) through dense sampling of each asteroid's light curve at each phase angle (equivalently, each epoch).

Using a large corpus of low-precision photometry from the MPC,⁹ Oszkiewicz et al. (2011, 2012) showed that a fitted parameter of one particular ϕ model correlates well with an asteroid's Sloan Digital Sky Survey (SDSS) visible color. While they were unable to correct for rotational variation (δ-term in Equation (1)), the Oszkiewicz et al. work nevertheless demonstrates a solid trend between ϕ and a compositional attribute (color).

These prior works motivate several defining aspects of this work's phase-function analysis:

• 1.  
  We fit multiple phase-function models to each light curve, both for compatibility with the literature and to explore how the fitted parameters are related.
• 2.  
  We simultaneously fit the rotational component with the phase-function part.
• 3.  
  We introduce a single colorimetric index for quantifying C-type versus S-type taxonomic classification, based on the compilation of several visible-band-color asteroid data sets (see Appendix), and examine the variation in phase-function parameters as a function of this color index.

The PTF is a synoptic survey designed primarily to discover extragalactic transients (Law et al. 2009; Rau et al. 2009). The PTF camera, mounted on Palomar Observatory's 1.2 m Oschin Schmidt Telescope, uses 11 CCDs (each 2K × 4K) to image 7.3 deg² of sky at a time at /pixel resolution. Most exposures (∼85%) use a Mould-R filter¹⁰ (hereafter "R"). The remaining broadband images acquired use a Gunn g-band filter. Nearly all broadband PTF images are 60 s integrations, regardless of filter. About 15% of nights (near full moon) are devoted to a narrowband (Hα) imaging survey of the full Northern Sky.

Science operations began in 2009 March, with a nominal one- to five-day cadence for supernova discovery and typical twice-per-night imaging of fields. Median seeing is 2'' with a limiting magnitude (for 5σ point-source detections), while dark conditions routinely yield (Law et al. 2010).

The PTF survey is ongoing and expected to continue through mid-2016. In January 2013 the PTF project formally entered a second phase called the intermediate PTF ("iPTF'; Kulkarni 2013). In this paper we simply use "PTF" to mean the entire survey, from 2009 through the present (2015). The iPTF program accommodates more varied "sub-surveys" as opposed to a predominantly extragalactic program, including variable star and solar system science. Images are still acquired with the same telescope/camera/filters with 60 s exposures and are processed by the same reduction pipeline.

Laher et al. (2014) describe the PTF data reduction and archiving pipelines, hosted at the Infrared Processing and Analysis Center (IPAC) at Caltech. Processing at IPAC includes bias and flat-field corrections, astrometric calibration against UCAC3 (Zacharias 2010), astrometric verification against Two Micron All Sky Survey (2MASS; Skrutskie et al. 2006), creation of source catalogs with Source Extractor (Bertin & Arnouts 1996), and production of reference images (stacks of ∼20–30 PTF images that reach ).

Ofek et al. (2012a, 2012b) describe the PTF survey's absolute photometric calibration method, which relies on source matching with SDSS DR7 (Abazajian et al. 2009), and thus requires at least partial overlap of PTF with SDSS each night. A separate, relative photometric calibration (based on light curves of non-variable field stars) also exists for PTF data and is described by Levitan et al. (2011) and in the Appendix of Ofek et al. (2011). In this work we utilize all R-band and g-band PTF data accumulated from the survey's start (2009 March) through 2014 July. The asteroid magnitudes reported in this work use relative photometric zero points when available (which as of this writing applies to ∼85% of PTF images) and absolute photometric zero points otherwise.

The PTF's robotic survey program and processing pipeline, as well as our data aggregation and analysis in this work, make use of many functions from the MATLAB package for astronomy and astrophysics (Ofek 2014).

Waszczak et al. (2013) used a custom spatial indexing algorithm to search the set of all PTF single-epoch transient detections (through 2012 July) for detections of all asteroids with orbits known as of 2012 August. That search procedure first generated uniformly spaced ephemerides for each asteroid using JPL's online service (HORIZONS; Giorgini et al. 1996). Each asteroid's ephemeris defines a 3D curve (two sky coordinates plus one time); the intersection of each curve with the 3D kd-tree of transient detections was then computed and positive detections within a 4'' matching radius saved.

In this work we use a modified version of the Waszczak et al. (2013) algorithm. The updates/changes are as follows.

First, in terms of content, we now search all PTF (R and g-band) data from 2009 March 01 through 2014 July 18 for all numbered asteroids as of 2014 July 12 (401,810 objects). We now exclude unnumbered objects as the positional uncertainty of these objects can be very large, and as they tend to be very faint, their light curves will not in general be of high quality.

Second, in place of a single-step matching of a 3D transient-detection kd-tree against 3D ephemeris curves, we now divide the search into two main steps. We first perform a 2D spatial matching that exploits the natural indexing of PTF exposures into tiles (i.e., the grid of evenly spaced boresights or "fields" on the sky). Each 2D ephemeris curve's intersection with the 2D PTF survey footprint is computed, the object's position is cubically interpolated to all epochs of exposures possibly containing the object, and the object's precisely computed positon is then compared to the precise image boundaries of candidate exposures. Matching of predicted positions against actual detections takes place subsequently as source catalogs are then loaded into memory (as needed and in parallel). This method is faster than the original Waszczak et al. (2013) method and enables separate logging of predicted and positive detections.

The results of the known-asteroid search, as well as the derived light-curve data (described later), are stored in a relational database, the size and contents of which are summarized in Table 1. Out of ∼18 million predicted single-epoch asteroid sightings (including predicted magnitudes as dim as , well below PTF's sensitivity), there were 8.8 million positive detections (within a 4'' radius). Of these, we define 4.3 million detections as "reliable" as they (1) lack any cataloged background sources within the 4'' radius, (2) have a calibrated magnitude uncertainty of less than 0.1 mag, and (3) lack any processing flags indicative of contamination. Figure 2 compares predicted, positive, and "reliable" detections; the middle and right panels of Figure 2 show that our definition of "reliable" seems to include a small fraction of likely bad observations (<1% contamination, note the vertical log scale), namely, those that have distance residuals greater than or magnitude residuals greater than ∼1 mag. Because these reliable detections are the subset of observations that we input into our light-curve-fitting model (Section 4), the fitting algorithm includes logic designed to remove isolated data points that have very large residuals, either with respect to the median light curve value or relative to their uncertainty.

Zoom In Zoom Out Reset image size

3.1.1. Intra-opposition Constraint

The most important parameter in the rotation component (the δ in Equation (1)) is the synodic spin period P, a constant that satisfies

where is the light-time-corrected observation timestamp, is the asteroid's geocentric distance, c is the speed of light, and n is any integer satisfying

where is the synodic orbital period,

where is the asteroid's sidereal orbital period and is its orbital semimajor axis (related by Kepler's third law). is the time elapsed between the asteroid's consecutive oppositions. Pursuant to this restriction, we constrain each δ solution using observations from within the same opposition—i.e., for most asteroids, within a 1.1–1.6 yr interval centered on the date of locally minimally observed α.

The intra-opposition restriction is important given that our data set (described in the next section) spans ∼5 yr. For an asteroid with a zero-inclination circular orbit and spin axis perpendicular to its orbital plane, we can relax Equation (3) to allow n to be any integer, in which case δ can be constrained using observations spanning many years. In general, however, Equation (2) must be modified to accommodate a varying viewing geometry with respect to the spin axis:

where F is some unknown periodic function satisfying , where m is any integer and is the sidereal orbital period. Provided that the amplitude of F is not large relative to that of δ, and provided that the spin vector is not changing with respect to the orbital plane (i.e., precessing¹¹ ) on a timescale comparable to , we are justified in assuming that Equation (2) (with the Equation (3) restriction) applies.

3.1.2. Second-order Fourier Series

Any δ satisfying Equation (2) can be approximated to arbitrary precision using a Fourier series. Harris et al. (2014) discuss why, from a geometric standpoint, the second harmonic tends to dominate an asteroid's fitted δ. As noted earlier (Section 1.1), most large asteroids approximately resemble triaxial prolate ellipsoids (e.g., Jacobi ellipsoids), having equatorial axis ratios of at most ∼3:1 (corresponding to a amplitude of ∼1.2 mag). For less extreme axis ratios (specifically, those producing a ∼0.4 mag or smaller second-harmonic amplitude), other harmonics related to shape or albedo asymmetries may contribute comparable coefficients to the Fourier approximation of δ.

The PTF survey program has—on a few rare occasions—conducted high-cadence (∼10 minutes spaced) observations of low ecliptic latitude fields. These runs produced a set of ∼1000 densely sampled main-belt asteroid rotation curves, which have already been analyzed and published (Polishook et al. 2012; Chang et al. 2014a). These high-cadence "pilot studies" are relevant to our present work in that they demonstrate (1) the quality of the PTF survey's photometric calibration for asteroids with unambiguously valid δ solutions and (2) the above-described prevalence of a dominant second harmonic in most of the objects sampled.

Following these pilot studies, we adopt a second-order Fourier series model:

where τ is the light-time-corrected epoch (see Equation (2)). In the pilot studies, most of the fitted δ solutions qualitatively resemble a simple sine or cosine function. Such a solution can be represented by either

• 1.  
  a first harmonic with period (with and ), or
• 2.  
  a second harmonic of period (with and ).

Given the prolate ellipsoid model, choice 2 is more realistic and hence preferred. However, again recognizing that other harmonics can have a non-negligible contribution, in fitting δ to our light-curve sample we allow the first-harmonic coefficients to be non-zero, but introduce logic into the fitting algorithm (see Section 4), which checks for double-period solutions satisfying certain criteria and iterates accordingly.

In this work we simultaneously fit each light curve's phase function ϕ along with its rotation curve δ (see Equation (1)). This approach is intermediate in complexity between some of the simpler, two-parameter (δ-neglecting) models that have been applied to very large data sets (e.g., Oszkiewicz et al. 2012; Williams 2012) and the more complex, shape plus pole-orientation models (Kaasalainen 2004; Cellino et al. 2009; Hanuš & Ďurech 2012) which can involve tens of parameters and require data spanning multiple oppositions.

Regarding the former class of models, we note that there is a formal statistical problem associated with neglecting δ when fitting ϕ. If modeling the observations M by , then the distribution of residuals is not Gaussian. Assuming that δ is a sinusoid with amplitude A, for observations M sampling the light curve at random times, the residual probability density function has a local minimum value at and maximum value near = ±A. Thus, p is bimodal and roughly bowl shaped—not at all Gaussian shaped. The uncertainty in ϕ produced by a standard minimization—which assumes Gaussian-distributed errors—is thus inaccurate. However, since p is symmetric about , for densely sampled data the fitted phase function ϕ remains unaffected by neglecting δ; in such a case the only effect is an underestimated uncertainty.

We obtain three separate fits for each light curve, each using a different phase function (ϕ) and allowing for unique solutions for H and δ in Equation (1). The three phase-function models are:

• 1.  
  the two-parameter model of Shevchenko (1997);
• 2.  
  the one-parameter G model (Bowell et al. 1989);
• 3.  
  the one-parameter G₁₂ model (Muinonen et al. 2010).

In this section we review and motivate the application of each of these ϕ models.

3.2.1. Two-parameter Shevchenko Model

Shevchenko (1997) introduced a phase function dependent on two parameters; in terms of Equation (1) the model is¹²

where β has units of mag/deg and C is the amplitude of the opposition surge (units of mag). This model was subsequently considered in-depth by Belskaya & Shevchenko (2000, hereafter B&S), who compiled the most complete (to date) set of high-precision, targeted phase-curve observations of main-belt asteroids from various data sets spanning several decades.

Though in practice Shevchenko's model is the least commonly used phase function out of the three we consider, it is by far the simplest to express mathematically, and it is the only model for ϕ whose parameters have linear dependence in Equation (1).

Furthermore, this model's parameters are the most straightforward to associate with physical asteroid properties. B&S highlighted a robust relationship between an asteroid's () phase-function parameters and its geometric albedo.¹³ As we later explore a similar relationship in the present work, we here review the basis of this observation.

The geometric albedo p_V is formally defined in terms of the phase function ϕ:

where is the (visible) bond albedo, defined as the total visible light energy reflected or scattered by the asteroid (in all directions) divided by the total visible light energy incident upon the asteroid (from the Sun). We also here define the phase integral q.

B&S showed that, in the range of β observed from S-type to C-type asteroids, β and C are empirically correlated, in a relation that we approximate here as

Using Equation (9) to substitute for C in Equation (7), inserting the result into Equation (8), and numerically evaluating the integral gives

B&S saw a negative correlation between p_V and β in the data,¹⁴ consistent with Equation (10) only if either is assumed constant among different asteroid types (not a reasonable assumption) or if negatively correlates with β, which B&S did not explicitly show.

The bond albedo can be thought of as an intrinsic, bulk-compositional characteristic of an asteroid's surface,¹⁵ much like an asteroid's color, whereas β and C relate (in part) to the textural, particulate, and macroscopic roughness of the asteroid's surface. B&S and other authors separately associate β with the shadow-hiding effect and C with the coherent backscatter effect. Both of these physical phenomena are understood from a theoretical standpoint (e.g., Helfenstein & Veverka 1989; Hapke 2012) to be functions of , with β negatively related to and C positively related. This is consistent with Equation (9) and renders Equation (10) consistent with B&S's noted p_V-versus-β correlation. Other properties such as particle size, particle geometry, and regolith porosity also have predicted (and laboratory-measured) contributions to the observed phase function (Hapke 2012 and references therein); these properties can conceivably vary independently of .

In short, our interpretation of the S-type and C-type asteroid data reviewed by B&S is that a compositional indicator () correlates with indicators of two independent phenomena (β and C) that contribute to how light scatters from an asteroid's surface. This statement intentionally makes no mention of p_V, since Equation (8) tells us that p_V by definition varies with β (in a non-obvious way) and with , the latter being a more basic compositional attribute.

As stated above, the phase function can be related to properties other than , such as regolith porosity. Many of these other properties in theory and experiment contribute to effects involving multiply scattered light and therefore do not alter the effect of shadow hiding (β-term in Equation (7)), which is dominated by singly scattered light (Hapke 2012). In contrast, the coherent backscatter effect (C-term) does involve multiply scattered light. B&S saw non-monotonic behavior in C as a function of p_V when including the rarer, high-p_V E-type asteroids in the same plot as C and S types. E types do conform, however, to the same negative monotonic trend in p_V-versus-β satisfied by the C and S types, consistent with the hypothesis that β is adequately expressed as a function of alone, yet E types have a lower-than-predicted C value based on extrapolation of Equation (9).

One possibility is that Equation (9) is not valid for all asteroids, but must be replaced by some unknown non-monotonic relationship, possibly because C depends non-monotonically on and/or has comparable dependence on other properties (e.g., porosity or grain size). Assuming that Equation (7) is a sufficiently general model for ϕ, and lacking knowledge of a good model for C, it follows that β and C should in practice always be fit separately. Another possibility is that Equation (7) is an incorrect or incomplete model; however, B&S described no instances wherein their model was unable to adequately fit the data for a particular asteroid or class of asteroids.

3.2.2. Lumme–Bowell G Model

The next phase-function model we consider is the Lumme–Bowell model (Bowell et al. 1989), also known as the (H,G) or IAU phase function:

Like Shevchenko's model, this model includes two terms (the basis functions and ) representing two physically distinct contributions to the observed ϕ. As detailed in Bowell et al. (1989), this model is semi-empirical in that it was derived from basic principles of radiative transfer theory with certain assumptions, and at various stages tailored to match existing laboratory and astronomical observations. That the two basis functions' coefficients are related to a single parameter G bears resemblance to the β-versus-C correlation described by Equation (9).

Marsden (1986) marked the IAU's adoption of this phase function as a standard model for predicting an asteroid's brightness. Since then this model has seen widespread application, and it is often used with the assumption G = 0.15 (e.g., in the ephemeris computation services offered by the MPC and JPL). Harris & Young (1988) present mean values of G for several of the major asteroid taxonomic classes (based on a sample of ∼80 asteroids), with G = 0.15 being an average between the C types () and the S types (). The G-model fails to accurately fit the rarer D types (which have linear phase curves) and E types (which have very sharp opposition spikes), whereas the Shevchenko model can properly accommodate these rarer types.

Use of the Lumme–Bowell ϕ in our light-curve model (Equation (1)) introduces a second nonlinear parameter (G) into the model, the period P being the other nonlinear parameter. This complicates the fitting algorithm somewhat, as described in Section 4.

3.2.3. Muinonen et al. G₁₂ Model

The third phase-function model we consider, introduced by Muinonen et al. (2010), bears resemblance to the G-model but includes a second free parameter and a third basis function:

As opposed to the analytic trigonometric basis functions of the G-model, here , , and (all functions of α alone) are defined in terms of cubic splines (see Muinonen et al. 2010 for the exact numerical definitions). Assuming that the coefficients G₁ and G₂ are constrained independently, these basis functions were designed to provide the most accurate fits to the phase functions of all major asteroid taxonomic types, including the rarer D types and E types.

For situations where fitting G₁ and G₂ separately is infeasible, Muinonen et al. (2010) specialized their above model to make it a function of a single parameter, G₁₂, which parameterizes G₁ and G₂ using piecewise functions:

In this work we use this single-parameter G₁₂ form of the Muinonen et al. model, making it analogous to the G-model in terms of implementation, including the complication associated with a nonlinear parameter.

3.2.4. Multi-parameter Hapke Model

Just as we commented on the more rigorous means of fitting a rotation curve via 3D shape modeling with multi-opposition data, for completeness we note that a more rigorous model (than the three presented above) exists for phase functions. Given better-sampled light curves and more computational power, future modeling of large photometric data sets would benefit from applying the more theoretically motivated model of Hapke (2012), an abbreviated form of which is

Here w is the single-scattering albedo (see footnote 15), of which r₀ is solely a function. The remaining factors all are functions of phase angle (α). Each opposition-surge term (B_S and B_C) has two free parameters (width and amplitude). K depends on the mean topographic roughness (a function of one free parameter); g is the single-scattering angular distribution function (typically includes one parameter); h is a function of α only; and is the phase function of an ideal Lambertian-scattering sphere (a simple function of α).

With its dependence, the Hapke model (Equation (14)) can conveniently eliminate both p_V and H from the modeling process. Inserting Equation (14) into Equation (1), and using the common relation¹⁶

where H is the absolute visual magnitude, D is the asteroid's effective diameter, and 1329 km is a constant (set by the arbitrarily defined magnitude of the Sun), produces a model with many physically meaningful parameters and free of both H and p_V.

For the particular case wherein we use Shevchenko's model (Equation (7)) for the phase function ϕ, the design matrix is

where the k-index has been omitted with the convention that k = 1 is the row of the above column vector, k = 2 is the second row, etc. Here and are the time and phase angle of the observation, f_j is the frequency, etc. Likewise, the linear-parameter matrix X in this case is

where H_j is the fitted absolute magnitude for the frequency, etc.

The general LLSq solution to Equation (16) is

where is the inverse of the data-covariance matrix B:

and is the parameter-covariance matrix, given by

where in the above definition we invert each of the matrices s_j, these being defined by

The elements of the parameter-covariance matrix S are the variances and covariances of the fitted parameters (as a function of frequency). The fit's residuals (as a function of frequency) are

and the fit's chi-squared (as a function of frequency) is

The frequency-dependent chi-squared is also known as the periodogram. Formally, the best-fit rotation frequency corresponds to the minimal value of , but this may differ from the preferred frequency solution if the light curve is contaminated by other systematic periodic signals, if the data suffer from underestimated measurement uncertainties, or if the best-fit frequency corresponds to a dominant first harmonic (as opposed to a preferred dominant second harmonic; see Section 3.1.2).

Figure 3 details our iterative light-curve-fitting algorithm's logic. Fitting commences as long as 20 or more "reliable" data points (see Section 2.2 and Figure 2) are associated with a light curve. Irrevocably bad data points are discarded in the first round of iterations; these include detections with 7σ or greater residuals from the initial solution. Examples of detections with such high residuals include contamination from background sources missing in the reference catalog, bad detector pixels that were not flagged by the pipeline, or spurious zero-point solutions.

Zoom In Zoom Out Reset image size

In the next stage of iterations, the fit's per degree of freedom is reduced to ∼1 (formally, it is reduced until it is less than 3; see Figure 3) by gradually inflating the observations' error bars through addition of a "cosmic error," so-named because it encompasses contamination from possible errors (in all the "cosmos"). In general, the cosmic error represents the same diverse contaminating phenomena responsible for the >7σ deviations seen in the initial iterations (see previous paragraph), just to a lesser extent.

Separately, this error bar inflation compensates for our model's inability to fit each asteroid's precise periodic structure using only two harmonic terms in the Fourier series. In the limit of infinite observations and sufficiently many Fourier terms, we would ideally expect our data's error bars to reflect true Gaussian variance. However, by truncating the series at two harmonics and using sufficiently precisely calibrated photometry, we are in effect choosing to sacrifice (downsample) some of our photometric precision to obtain a formally better fit at the coarser resolution limit of the model.

To illustrate use of the cosmic error, consider the example of an eclipsing binary light curve, i.e., a rotation curve that is effectively sinusoidal except for a small interval around the phase of minimum flux, when it dips to a lower-than-predicted brightness. Observations acquired during such eclipses will have systematic negative deviations greater in absolute value than would be explained by Gaussian variance alone. Increasing the error bars of these observations will decrease the fits' without altering the value of the fitted frequency. The fitted parameters' uncertainties (for both frequency and the linear parameters) are accordingly inflated as a penalty, and the fitted amplitude will be underestimated. As detailed in Figure 3, the initial cosmic error used is 0.002 mag, and each iteration is multiplied by a factor of 1.5 until the is sufficiently low. If the cosmic error exceeds 0.1 mag, the fitting is aborted. If the (per degree of freedom) drops below 3 while the cosmic error is still below 0.1 mag, the fitting process concludes "successfully" (see Figure 3).

Concurrently, each iteration includes a test for the presence of double peaks in the folded rotation curve (only if the fitted amplitude is at least 0.1 mag). In particular, if there exist two maxima and two minima in the folded light curve, we demand that the ratio of these peaks be greater than 0.2. Such a solution is preferred (see Section 3.1.2) given our ellipsoidal shape assumption, as described by Harris et al. (2014).

Denote as the frequency yielding the absolute minimum per degree of freedom value, denoted (after the cosmic error has been tuned). If the folded light curve is single-peaked (or has only a relatively small secondary peak), then another deep minimum usually exists at the harmonic frequency , the local minimum value of which we denote . For cases wherein , where is the inverse of the cumulative distribution function for N free parameters evaluated at p, then we instead choose rather than . The 1σ uncertainty interval for the best-fit frequency is then found by computing the upper and lower intersections between and the periodogram in the vicinity of . Note that we used n = 7 free parameters in this case, i.e., the number of elements of X_j (Equation (18)).

Modeling the phase function ϕ with either the G or G₁₂ model (Equations (11) and (12)) introduces a second nonlinear parameter (after the frequency f), and so we must modify the equations of the previous section accordingly. We sample evenly spaced phase-function parameter values. In particular, for G we test the interval in steps of , and for G₁₂ we test the interval in steps of .

Our approach is to modify the left-hand side of Equation (16) by defining a new matrix that contains all possible phase-function-corrected observed magnitudes:

where, e.g., for the case of the G-model (Equation (11)),

The linear-parameter-solution array X now has an extra index q, reflecting the fact that we are now solving for each linear parameter as a function of the two nonlinear parameters. The design matrix has the same number of indices as before (but fewer rows):

while the linear-parameter matrix X is now

The appeal in adopting the above approach is that the general solution is only slightly modified:

where the only differences between Equations (19) and (29) are the q indices appended to X and m (and the latter being redefined as m').

The fit's residuals R are now a function of frequency and phase-function parameter:

as is the fit's chi-squared:

As a function of any of the linear parameters, the fit's varies precisely quadratically, whereas as a function of frequency it has an intricate spectral structure with many local minima. As a function of a nonlinear phase parameter (G or G₁₂), the tends to have a single minimum (on the range we evaluate): in this sense G and G₁₂ are more similar to the linear parameters than they are to frequency. However, the generally asymmetric shape of the phase parameter' dependence necessitates its grid-based numerical treatment—particularly to ensure accurate estimation of the phase parameter's uncertainty.

The two-dimensional surface given by Equation (31), which is defined on an grid, can be reduced to a one-dimensional function by choosing, for each frequency index j, the phase-parameter index q that minimizes the . The result is a one-dimensional periodogram, as in Equation (24). Once the fitted frequency is identified, we compute the uncertainty in the fitted f by the method described in the previous section using the inv--cdf() function. We then likewise numerically compute the uncertainty in the phase parameter by again collapsing to a one-dimensional vector, this time as a function of the phase parameter with the frequency fixed at the fitted value (j-index), and use the inv--cdf() function to estimate the uncertainty in the phase parameter.

As noted in Table 1, a total of 587,466 light curves exist in PTF, where each light curve by definition consists of all reliable observations of a unique asteroid observed in a single opposition in a single photometric band. Of these, only ∼10% (59,072 light curves) have at least 20 observations and therefore qualified for fitting with our algorithm. A total of 54,296 light curves actually produced a fit—the remaining ∼5000 light curves failed to produce a fit either because some observations were discarded and the total fell below 20 data points or because the fitted cosmic error grew to exceed 0.1 mag.

Figure 4 shows several examples of light curves fitted with the algorithm described in this section. In the third column (column C) of Figure 4, we show the periodograms of each light curve. Note that although the periodogram's horizontal axes are labeled with the period (for easier interpretation), the chi-squared (per degree of freedom) values are actually plotted linearly with respect to frequency. This is because, as described earlier, our sampling is uniform with respect to frequency, and the harmonics are more easily seen with constant frequency spacing. Column (D) shows the data sampling in rotational phase versus solar phase angle, a useful plot to ensure that there is no obvious correlation between the two (which could lead to an erroneous fit, e.g., for long periods, large amplitudes, and/or few data points).

Zoom In Zoom Out Reset image size

Each iteration in the fitting of each asteroid light curve involves evaluating the arrays and tensor-products in either Equation (19) or (29). This includes inverting the data-covariance matrix B (Equation (20)) and inverting the matrices s_j (Equation (22)). The arrays L, m', X, and R can have a relatively large number of elements, making them and their relevant products potentially taxing with respect to computational memory.

Our particular implementation of this algorithm leverages the efficient array-manipulation capabilities of MATLAB, especially its ability to perform fast matrix multiplication and matrix inversion utilizing BLAS calls¹⁷ and OpenMP multithreaded C loop code.¹⁸ Given typical numbers of observations and frequency sampling, each of our light-curve fits (including the multiple iterations) takes on average several tens of seconds to run on an eight-core machine (multithreading enabled) and typically consumes less than ∼4 GB of memory using single-precision computation.

In a tar.gz file we provide our custom MATLAB function used for fitting the G-parameter version of the light-curve model (asteroid_lc_fit_G.m). Analogous versions exist for the Shevchenko and G₁₂ models. This function takes as input an asteroid's apparent magnitudes, magnitude uncertainties, observed epochs, phase angles, and geocentric and heliocentric distances. Its outputs include the linear-parameter-solution array (Equation (28)), residuals (Equation (30)), chi-squared array (Equation (31)), and additional information about each light-curve solution such as the amplitude and peak ratios.

A total of 927 (∼2%) of our fitted light curves belong to 805 unique asteroids having a previously measured period listed in the LCDB of Warner et al. (2009). This includes only asteroids having a quality code of 3 (highest quality) in the LCDB.

Figure 5 shows that the distribution of relative errors on our fitted frequencies is bimodal, with the left mode corresponding to periods having better than ∼3% agreement with the reference period, and the right mode corresponding to periods in disagreement with the reference period. These disagreeing fits include light curves that differ from the reference value by a harmonic (half = relative error 0.5, double = relative error 1.0), as well as frequencies that do not differ by a factor of two or any integer multiple. About 1/3 of the light curves in Figure 5 fall into the right mode and are thus considered disagreeing fits.

Zoom In Zoom Out Reset image size

Figure 6 shows some examples of these disagreeing fits. Row 1 shows an apparent low-amplitude rotator, whose fitted period of 15.7 hr differs from the reference value of 9.7 hr. Row 2 is an object whose periodogram contains a great deal of noise, divided into two broad forests of frequency minima. The left forest appears to have been selected by our fitting algorithm, while the right forest seems associated with the true period of ∼2.7 hr. Row 3 contains an object whose 12% relative frequency error exceeds the 3%-accuracy threshold we have defined, and so despite appearing to be a good fit it is formally categorized as inaccurate. Row 4 also looks like a reasonable fit at 6.4 hr, but disagrees with the reference period of 11.0 hr (though the latter does have a perceptible local minimum in the periodogram). Finally, Row 5 includes a likely example of the algorithm fitting noise in the photometry of a faint asteroid.

Zoom In Zoom Out Reset image size

In Figures 7 and 8 (top and middle rows) we detail the distribution of the accurately recovered-period and inaccurately recovered-period subgroups in terms of eight different light-curve parameters. Some basic observations from these histograms are:

• 1.  
  fitted periods are far less reliable if longer than ∼1 day or shorter than ∼2.7 hr;
• 2.  
  fitted amplitudes of less than 0.1 mag correspond to the least reliably fit periods;
• 3.  
  light curves consisting of observations dimmer than ∼18.5 mag are much less reliable than brighter light curves (though they are also far less numerous in the known-period sample);
• 4.  
  fit (per degree of freedom) values of less than ∼1.7 correlate with less reliable periods (though they are also far less numerous in the known-period sample); note that, in the fitting process, growth of the cosmic error term ceased once the (per degree of freedom) fell below 3 (see Figure 3);
• 5.  
  the number of observations in a light curve is not directly correlated to the reliability of the fitted period;
• 6.  
  the ratio of the folded light curve's two peaks, the signal-to-noise ratio of the periodogram's chosen minimum, and the uncertainty in the absolute magnitude parameter are all strong indicators of the reliability of the fitted period.

Zoom In Zoom Out Reset image size

The above comments reflect consideration of the one-dimensional distributions in Figures 7 and 8; however, we can easily imagine that there are correlations in more dimensions not evident from these plots alone. An obvious example would be the two-dimensional distribution in amplitude versus median magnitude: reliability is presumably greater for bright asteroids having amplitudes mag than it is for dim asteroids having amplitudes mag. Period versus amplitude is also likely an insightful distribution (and was considered, for example, by Masiero et al. 2009). The number of observations possibly would correlate with reliability if we were to restrict another parameter or parameters to some specific interval.

Zoom In Zoom Out Reset image size

Rather than manually examining the period-fitting reliability as a function of all possible multi-dimensional combinations of the eight light-curve parameters detailed in Figures 7 and 8, we can take a more general approach of considering the reliability to be a single function defined on the multi-dimensional parameter space in which all of the light curves reside. We hypothesize that accurately fit light curves and inaccurately fit light curves occupy distinct regions in this multi-dimensional volume. As these volumes can overlap to some extent, we can at least estimate the probability that a light curve with that particular vector of parameters corresponds to an accurately recovered (or inaccurately recovered) period when obtained by the fitting algorithm of Section 4.

There are two general ways of accomplishing this goal. One way is to produce a large number of synthetic light curves filling out the multidimensional light-curve-parameter space, subject these synthetic light curves to our fitting algorithm, and thereby map out, e.g., by binning and interpolation, the fit reliability throughout the multi-dimensional volume. This method requires us to accurately simulate all sorts of varying sampling cadence as well as measurement uncertainties, including contributions from both systematics and noise, and it requires significant extra computing time to actually subject the synthetic data to our fitting procedure. The second method—the approach we take in this work—uses a ground-truth sample (the known-period light curves already described in this section) to train a machine classifier to discriminate reliable versus unreliable fits within the multi-dimensional light-curve-parameter space.

We adopt a supervised ensemble-method approach for classification, originally popularized by Breiman et al. (1984), specifically the random forest (RF) method (Breiman 2001). RF classification has extensive and diverse applications in many fields (e.g., economics, bioinformatics, sociology). Within astronomy in particular RF classification is one of the more widely employed methods of machine learning, though many alternatives exist. For example, Masci et al. (2014) use the RF method for variable-star light-curve classification, while others have approached this problem via the use of, e.g., support vector machines (Woźniak et al. 2004), Kohonen self-organizing maps (Brett et al. 2004), Bayesian networks and mixture-models (Mahabal et al. 2008), principle component analysis (Deb & Singh 2009), multivariate Bayesian and Gaussian mixture models (Blomme et al. 2011), and thick-pen transform methods (Park et al. 2013).

For general descriptions of RF training and classification, we refer the reader to Breiman (2001), Breiman & Cutler (2004), and the many references cited by Masci et al. (2014). Our use of an RF classifier is particularly motivated by its already-proven application to the discovery and classification of astrophysical transients in the same PTF survey data (Bloom et al. 2012), as well as streaking near-Earth asteroid discovery in PTF data (A. Waszczak et al. 2015, in preparation).

Machine-learning application generally consists of three stages: training, cross-validation, and classification. In the training stage of building a machine classifier, the multi-dimensional parameter space is hierarchically divided into subspaces called nodes; these nodes collectively comprise a decision tree. The smallest node—also known as a leaf—is simply an individual data point (in our case, a single light curve). Given a set of leaves with class labels, one can build an ensemble of trees (called a forest), each tree representing a unique partitioning of the feature space, wherein the nodes are split with respect to different randomly chosen subsets of the parameter list. Each node splitting attempts to maximize the separation of classes between the sub-nodes. Serving as a model, in the subsequent classification stage the forest allows one to assign a probability that a given vector of features belongs to a given class. During cross-validation (an essential early stage in this process), the training and classification steps are repeated many times, each time using different subsamples (of labeled data) as the training data and testing data. Cross-validation evaluates the classifier's performance and ensures that it is not overfitting the training data.

For our light curves, we are interested in a binary classification, i.e., whether the fitted period is accurate ("real") or inaccurate ("bogus"). Bloom et al. (2012) coined the term realBogus to describe this binary classification probability in the context of extragalactic transient identification. In the present work we are essentially adapting Bloom et al.'s realBogus concept to the problem of light-curve-period reliability assessment.

We employ a MATLAB-based Random Forest classifier,¹⁹ which is a port of the original RF software (originally implemented in R). This software includes two main functions, which perform the training and classification steps separately.

5.2.1. Classifier Training

Our training data consist of the known-period light curves (see the previous section) belonging to the two classes under consideration: 618 light curves having accurately fit rotation periods and 309 light curves having inaccurately fit periods. Membership in one class versus the other depends on our arguably arbitrary 3% relative accuracy threshold, though we claim that the clearly bimodal shape of the distribution in Figure 5 justifies this 3% criterion. We note also that the classifier ultimately only provides a probability that a given light curve belongs to one class or the other, so that objects very near to the 3% cutoff may conceivably correspond to classification probabilities close to 0.5.

An important point is that the "ground-truth" reference periods we have taken from the database of Warner et al. (2009) may include some number of inaccurate periods. Such periods may be the product of erroneous fitting on the part of any one of its many different contributors, each of whom may employ a different fitting procedure and/or adhere to different confidence criteria. For the sake of this work, however, we consider all quality code 3 periods to be accurate—any actual inaccuracy will contribute to decreased classifier performance.

Besides ground-truth periods that are simply inaccurate, we also, in principle, risk contamination from reference periods that are no longer accurate. We assume that the majority of asteroids' periods are not changing with time, at least not at levels measureable with our data. For instance, direct measurement of the YORP mechanism in at least one asteroid (Lowry et al. 2007) reveals a relative rotation period change of several parts per million over several years. Any measureable period changes would likely be due to recent collisional events. The case of asteroid 596 Scheila (Bodewitts et al. 2011) demonstrates that detectable collisional events among main-belt asteroids do occur on a relatively regular basis, though even this robustly detected collision imparted no measurable change in the asteroid's spin rate (Shevchenko et al. 2013).

Although Figures 7 and 8 detail the period-fitting reliability as a function of only eight light-curve parameters, we construct our classifier using 12 additional parameters, for a total of 20 light-curve parameters. In the context of machine learning these parameters are known as features. The 20 features we use were chosen on the basis of their availability (most are output directly by the fitting process and do not require additional computation), as well as their actual importance (as computed during the cross-validation tests described in the next section).

Figure 9 visualizes the two-dimensional correlation coefficients for all possible pairs of the 20 light-curve features. Overall, the correlation structure of the training sample qualitatively resembles that of the full data set, implying that the training set fairly well represents the overall data set in terms of its feature-space structure. On the other hand, the distributions (e.g., median value, range of values) of individual features in the training set do not necessarily match the distributions in the full data set. An obvious example is that the full data set contains far more faint asteroids than does the training sample, even though in both cases the median magnitude (medMag) is positively correlated with quantities like rmsFit (due to Poisson noise) and hMagRef (since larger asteroids tend to be brighter).

Our 20 light-curve features are listed in Table 2, in order of decreasing importance. Most of these quantities we have discussed already in previous sections in the context of our model and fitting procedure. The list also includes two features characterizing the magnitude distribution of the folded light curve: (1) Stetson's K-index, a measure of the kurtosis borrowed from variable-star-light curve analysis (Stetson 1996), and (2) a "cusp index" that quantifies the extent to which the dimmest 10% of the data points in the folded light curve deviate from the best fit relative to the other 90% of the data points. We designed the cusp index to potentially identify eclipsing systems that are poorly fit by the two-term Fourier approximation but nonetheless may have accurately fit periods (examples of light curves with such cusp-like minima appear in Figure 10). Eclipsing binaries would be most properly treated with a different model entirely, as would tumbling asteroids (which we also did not systematically try to identify in the data, and probably lack reliable light curve solutions when subjected to this work's algorithm).

Table 2.  Summary of the 20 Light-curve Parameters (Features) Used by Our Period-Quality Classifier

Feature	Importance (%)	Description
peakRatio	11.1	Ratio of the fitted light curve's two peaks (). Zero if only one peak, one if exactly the same height.
amplitude	10.2	Fitted amplitude of the folded light curve. Equivalent to the height () of the larger of the two peaks.
periodFit	8.6	Rotation period value obtained using this work's data and fitting algorithm.
freqSNR	8.4	Signal-to-noise ratio of the fitted (minimum) frequency in periodogram =
hMagErr	5.8	Uncertainty in the fitted H-magnitude (i.e., error in the fitted absolute magnitude)
a12Coeff	4.3	Fourier coefficient A12
a22Coeff	4.2	Fourier coefficient A22
numObsFit	4.1	Number of observations in the final fitted light curve, after discarding any bad observations
medMag	4.1	Median calibrated magnitude (in the photometric band specific to the light curve, either R or g)
chisq	4.1	Reduced chi-squared of the fit, i.e., per degree of freedom)
a21Coeff	4.0	Fourier coefficient A21
a11Coeff	3.8	Fourier coefficient A11
rmsFit	3.8	Root-mean-squared residual of the fit
hMagRef	3.8	Reference H-magnitude (i.e., absolute magnitude of the asteroid in V band as listed by the MPC)
kIndex	3.7	Stetson's K-index (a measure of kurtosis in the magnitude distribution of a folded light curve, introduced by Stetson 1996)
freqResol	3.7	Resolution of the periodogram: , where is the time between the first and last observations in the light curve
hMagResid	3.7	Difference between the reference absolute magnitude (hMagRef) and the fitted H-magnitude
cuspIndex	3.6	"Cusp index": Median squared residual of the dimmest 10% points divided by the median squared residual of all other points
numObsRem	2.9	Number of observations removed during the fitting process (due to >7σ residuals with respect to preliminary fits)
cosmicErr	2.1	Final "cosmic error" value at end of fitting process ( mag in all cases)

Note. See text for a discussion of the cross-validation-derived importance value (Section 5.2.2).

Download table as:  ASCIITypeset image

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

5.2.2. Classifier Cross-validation

To ascertain the trained classifier's capabilities, and to ensure that the classifier is not overfitting the training data, we perform a series of 1000 cross-validation trials. In each trial we split each class (accurate fits and inaccurate fits) into a training subsample (a randomly chosen²⁰ 80% of the class) and a test subsample (the remaining 20% of the class). We then train a classifier using the combined training subsamples and subsequently employ the classifier on the combined test subsamples. In each of the trials, the classifier outputs a classification probability (score) for each object in the test sample, and we track the TPR (fraction of accurate period fits that are correctly classified above some threshold probability) as a function of the FPR (fraction of inaccurate period fits that are incorrectly classified above said threshold probability). See Figure 11 for a summary of these terms.

Zoom In Zoom Out Reset image size

The results of the cross-validation are shown in Figure 12. By tuning the minimum classification probability used to threshold the classifier's output, one effectively moves along the hyperbola-shaped locus of points in TPR-versus-FPR space seen in the plot. Several points have labels (p = ...) indicating the corresponding threshold probability (adjacent points being separated by ). The error bars in Figure 12 represent the standard deviation of the location of each point over all 1000 trials, while the point centers are the average locations.

Zoom In Zoom Out Reset image size

A classification threshold of is conventionally used when quoting single FPRs and TPRs. In our case, this gives FPR with TPR . The contamination of positively classified light curves in the cross-validation trials depends also on the actual class ratios in the sample being classified. In particular, since of our known-period light curves are inaccurate fits (Figure 5), it follows that among all light curves the classifier labels as accurate fits, the contaminated fraction is . If instead of using the classifier we just randomly labeled some fraction of the light curves as accurate and the rest as inaccurate, the resulting contamination would be (i.e., worse than the afforded by the classifier, as expected).

Several parameters can be adjusted or tuned when training an RF classifier. First is the number of decision trees generated during the training stage. Classification accuracy typically increases with the number of trees and eventually plateaus. Most applications employ hundreds to thousands of trees; we here use 1000 trees. Another tunable parameter is the number of randomly selected features (out of the 20 total here considered) with respect to which nodes are split in building the decision trees. Breiman (2001) recommends using the square root of the number of features. We ran the cross-validation for all possible numbers of features with respect to which the nodes can be split (i.e., all numbers between 1 and 20). The results are in Figure 13. We chose n = 4 as the number of features to split, both because the classifier's performance plateaus after that value and because it follow's the recommendation of Breiman (2001) () features.

Zoom In Zoom Out Reset image size

Other parameters that can be tweaked are the maximum depth of a tree, the minimum number of samples per leaf, the minimum number of samples used in a split, and the maximum number of leaf nodes. We do not constrain any of these parameters, meaning we allow trees of any depth, with any number of leaf nodes, leaf nodes consisting of a single sample, and splits based on the minimum of two samples. We note that as a result our model optimization is not comprehensive and it is possible that a better classifier could be trained. However, the relatively small training sample size here is likely the limiting factor; additional data are necessary to substantially improve the classifier performance.

In the bottom rows of Figures 7 and 8, we detail the dependence of the TPR and FPR on various light-curve parameters. Averaging (marginalizing) over any of the x-axis quantities in these bottom-row plots (while also weighting each bin by the number of light curves it contains; see the top row of plots in Figures 7 and 8) produces precisely the TPR and FPR values of the p = 0.5 data point in Figure 12.

In addition to the TPR and FPR estimates, cross-validation allows us to quantify the relative importance of the features by computing the average depth in the trees at which a split was performed with respect to each feature. Those features with respect to which the training sample is consistently divided early in the building of each tree are deemed more important (i.e., more discriminating) than those features that are split later, as the tree-building process tries to maximize the separation of the classes as early as possible by splitting features in an optimal sequence. Both Table 2 and Figure 9 list the features in order of importance.

Note that we had manually guessed several of the most important features—namely, peakRatio, freqSNR, and hMagErr—prior to any machine-learning work via inspection of the plots in Figure 8. The numerical importance values thus agree with these initial observations and also quantify the significance of features that would be difficult to ascertain manually. For instance, numObsFit appears (in Figure 8) not to be related to the fitting accuracy, while medMag (Figure 7) does appear related to accuracy (fainter light curves being less accurate), yet these two features evidently have equal importance in the classification process (see Table 2). Figure 9 indicates that numObsFit and medMag have quite different correlation relationships with respect to more important features. Hence, it would not be surprising if their one-dimensional distributions (in Figures 7 and 8) bear no resemblance to the multi-dimensional distributions on which the decision trees are defined and in which these two parameters apparently carry comparable weight.

5.2.3. Machine-vetted Light Curves

In addition to machine-based vetting, we manually inspected all 54,296 of the light curves that were output by our fitting process. A human screener first studies the ground-truth known-period examples (Section 5.1) in an effort to learn to distinguish between accurate and inaccurate fits. Only the G₁₂ fit is considered (as was the case with the automated classifier), and for each light curve the screener inspects precisely the amount of information included, for example, in Figures 4, 6, and 10 of this paper. Specifically, for each light curve the screener views a row of four plots: (1) the rotation-corrected phase curve, (2) the phase-function-corrected folded rotation curve, (3) the periodogram, i.e., the reduced plotted linearly against frequency (labeled however with the corresponding period), and (4) the rotational-phase versus phase-angle plot. A single screener is presented with these plots through a plain-formatted webpage, allowing for efficient scrolling through the light curves and rapid recording of either a "reliable" or "unreliable" rating for each fitted period. In addition, all light curves in the known-period sample were reinserted into the screening list, with their reference periods removed. These were thus blindly assessed by the screener, independent of their formal (3% accuracy) classification status.

The black lines in Figure 14 plot the results of the manual screening, in which a total of 10,059 light curves (19% of the total set) were deemed "reliable." With respect to the machine-approved sample, the human-rated sample is in all cases between roughly a factor of ∼1 and 2 smaller in each bin relative to the features examined in Figure 14. In general, the shapes of the machine-approved and human-approved distributions match fairly closely. Figure 15 shows examples of light curves for which the machine- and human-based classifiers differed in their rating (we focus on very short and very long fitted periods in Figure 15, but many examples exist for intermediate periods as well).

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

A total of 654 unique asteroids have more than one PTF light curve whose fitted period was labeled as reliable by the vetting process described in the previous sections. These 654 asteroids collectively have 1413 fits (so that the average multiplicity is ∼2.2 fits per asteroid) and include objects either observed in multiple oppositions and/or in both filters during one or more oppositions. Figure 16 plots the distribution of the relative error in the fitted frequencies of all such multiply fit asteroids, this error being defined as the range of the asteroid's fitted frequencies divided by the geometric mean of its fitted frequencies. Just as in Figure 5 (when we compared to literature-referenced frequencies), we see a prominent mode in the histogram peaking at ∼0.1% relative error, with some excess for errors greater than ∼3% error. There are 63 asteroids in particular with relative errors greater than 3%; of these, only four asteroids have more than two fits. If we assume that, in the remaining 59 pairs of disagreeing periods, one of the periods is correct, then the contamination fraction of light curves based on the sample of multiply fit asteroids is %.

Zoom In Zoom Out Reset image size

The preliminary demographic analyses that follow do not take into account fully de-biased distributions of, e.g., spin rates, amplitudes, or phase-function parameters. The TPRs and FPRs given in the bottom row of plots in Figures 7 and 8 (also, the blue and violet lines in Figure 14) constitute some of the necessary ingredients for producing a fully de-biased data set; however, in this work we do not attempt to compute the de-biased distributions.

In Figure 17 we reproduce several of the plots appearing in Pravec et al. (2002 and references therein), using this work's much larger data set (characterized by at least an order of magnitude larger sample of small objects). Both spin rate and amplitude are examined for the 4040 objects having diameter data from infrared surveys. Unlike Pravec et al. (2002), we are not able to individually plot each light curve's data (the ∼4000 points would make the plot difficult to render, as well as difficult to read); hence, we plot these (and other relationships later in this section) using two-dimensional histograms where the intensity of each pixel corresponds to the number of objects in that bin (darker means more, with linear scaling). Additionally, 2D histograms for which the diameter is plotted on the horizontal axis have their pixel values column-normalized, i.e., all pixels in each column of the histogram sum to the same value. This facilitates the visual interpretation of period and amplitude variation with diameter, as the left-hand side (small-diameter end) of the plots would otherwise saturate the plot.

Zoom In Zoom Out Reset image size

Following Pravec et al. (2002), we include the geometric mean rotation frequency as computed from a running bin centered on each object. The half-width of the bin centered on each object is either 250 (data points) or the object's distance from the top or bottom of the sorted diameter list, whichever is smallest. This ensures that the geometric mean is not contaminated at the edges of the plot by the interior values, though it also means that more noise exists in these edge statistics. The geometric mean is the more intuitive statistic for the rotation period as compared to the arithmetic mean, since the rotation periods tend to span several orders of magnitude. In addition to the geometric mean, we plot the 16th and 84th percentile values from each running bin.

The basic observed trend regarding rotation rate is that smaller-diameter asteroids rotate faster on average. A slight increase in the rotation rate also appears for objects larger than ∼80 km. Binning the data into a coarser set of three diameter bins and normalizing each object's spin rate by the local geometric-mean rate, we see a progression from a near-Maxwellian distribution to a progressively non-Maxwellian distribution for smaller objects. The rotation rates of a collisionally equilibrated population of rotating particles are known to approach that of a Maxwellian distribution (e.g., Salo 1987), which for a population of N objects as a function of rotation frequency f is

where is the number of objects in the interval and is the peak frequency (i.e., the frequency corresponding to the distribution's maximum).

One way of testing how well a Maxwellian actually fits the data is the two-sided Kolmogorov–Smirnov test (Massey 1951). This test compares an empirical distribution to a reference distribution (e.g., Gaussian, Maxwellian, or another empirical sample) via a bootstrap method. In particular, it computes a statistic quantifying the extent to which the cumulative distribution function differs in the two distributions being compared. In our case, we use Equation (32) to simulate a large sample (10⁵) randomly drawn from an ideal Maxwellian distribution and compare this simulated sample against the 99-asteroid sample (of km) frequencies. Interestingly, this test indicates that our 99 large-asteroid normalized frequencies differ from a Maxwellian at nearly the 10σ confidence level, though this could be due in part to the lack of a proper de-biasing of the distribution (see Section 6.1)

All of these trends—including the qualitative resemblance of a Maxwellian but its formal disagreement—were noted previously by Pravec et al. (2002). At the time their km size bin contained data on only 231 objects, as opposed to our sample of 2,844 asteroids with km. Conversely, our km bin contains only 99 objects as compared to the ∼400 large asteroids they took into consideration in comparing to a Maxwellian.

Steinberg & Sari (2015) recently described how collisional evolution of large asteroids should actually lead to a Lévy distribution, which has a significantly longer tail than a Maxwellian distribution having the same peak. They compared their theory to spin rates of km asteroids from the LCDB and found in general that the Lévy distribution fails to fit the spin distribution of large asteroids, suggesting that there may be a significant primordial component to the spin distribution. Potential primordial contributions to the angular momentum of asteroids were explored by Harris & Burns (1979) and later authors.

Our amplitude distribution contains an obvious observational bias (see Section 6.1) in that amplitudes less than ∼0.1–0.2 mag are generally ill fit by our modeling procedure (see Figure 7) and thus significantly underrepresented in our sample of reliable light curves considered here. Nonetheless, we see a clear trend of smaller asteroids exhibiting larger rotational amplitudes, consistent with the idea that larger bodies have sufficient surface gravity to redistribute any loose mass to a more spherical shape.

As we have done for the normalized frequency distribution, we plot diameter-binned normalized amplitudes against a Maxwellian distribution, this time merely to guide the eye as opposed to validating any hypothetical physical interpretation. The fact that the normalized amplitude distributions do not deviate too drastically from the Maxwellian shape at smaller diameters indicates that the spread in the amplitude distribution is proportional to its mean value, a basic property of the Maxwellian distribution, hence the good agreement. Carbognani (2010) provides a recent analysis of asteroid rotation amplitudes and highlighted a similar increase in both the amplitude's mean and spread with decreasing diameter.

Panel (B) of Figure 17 shows the distribution in period-versus-amplitude space, in which we can plot all 9033 light curves, including those lacking a diameter estimate. Contours representing the maximal spin rate of a body held together solely by self-gravity of certain uniform densities are overplotted. Our data as a whole do not appear to populate the region beyond the ∼2 g cm⁻³ contour. Later in this section we will re-examine this behavior separately for the two major taxonomic classes.

We consider any of the 54,296 fitted PTF light curves to have a reliably fit phase function if both of the following conditions are satisfied:

• 1.  
  The light curve is one of the 9033 having a reliable period fit, or its fitted amplitude (for the G₁₂ model) is less than 0.1 mag (the latter is true for 1939 light curves, only 39 of which have reliable periods).
• 2.  
  The light curve is fit using data from at least five phase-angle bins of width . These five bins need not be contiguous, and they need not include phase angles in the region where opposition surges are typically measured (i.e., ).

The above two criteria are met by 3902 out of the 54,296 PTF light curves. Of these, 1648 have an infrared-based diameter available, 651 have a color index available, and 361 have both a diameter and color index.

Figure 18 details the distributions of the fitted phase parameters G₁₂, G, β, and C against the color index, bond albedo, and in 1D histograms with color-based taxonomic subsets. Though the phase parameters are all correlated with color index and with bond albedo, none of the 1D phase-parameter distributions (right column of plots) exhibit bimodality alone, whereas the bond albedo (bottom right plot) does show significant bimodality. The red and blue histograms consist of all asteroids having color metric either less than 0.25 (C types) or greater than 0.75 (S types). The G and phase parameters are only plotted for those light curves that also have a G₁₂ solution. Not every light curve produced a solution for all three of the phase-function models; hence, the sample sizes for the G and models include a slightly reduced number of data points.

Zoom In Zoom Out Reset image size

We reiterate our statement from Section 3.2.1 that the bond albedo is a more fundamental (i.e., intensive rather than extensive) property than is the geometric albedo p_V, hence our focus on here. The bond albedo is computed using Equation (8) together with Equation (15) and makes use of our PTF-derived absolute magnitudes—H from the G₁₂ fit in particular—as well as the phase integral q of Equation (8), also computed directly from the G₁₂ fit for ϕ. In particular,

6.3.1. Taxonomy from Light-curve Data

We use the distribution of bond albedo versus G₁₂ to define another taxonomic metric analogous to the color index. In particular, we apply the same clustering analysis to this distribution as we did for the seven 2D color distributions in the Appendix. This procedure assigns to every object in the -versus-G₁₂ diagram a probability of membership in each of two clusters (color-coded blue and orange in Figure 19). The cluster centers are fit by the algorithm, and the output class probability of a given data point relates to its distance from these cluster centers. Probabilities near 0 represent likely C-type class membership, while probabilities near 1 represent likely S-type membership. We refer to this new metric as the photometric index; it complements the color index as another proxy for taxonomy. There are 361 asteroids with both a photometric index and color index available (Figure 19, right plot); the two indices are clearly correlated (, significance). Note that asteroids only have a defined photometric index if they have an infrared-derived diameter available, so that is defined.

Zoom In Zoom Out Reset image size

6.3.2. Wavelength Dependence

Observational evidence for the reddening of asteroid colors with increasing phase angle is discussed by Sanchez et al. (2012 and references therein). Color variation with phase angle can be equivalently stated as variation of the phase function with wavelength. Asteroids that have PTF light curves in both of the survey's filters (R and g band) allow us to investigate this phenomenon. We note, however, that Sanchez et al. (2012) describe phase reddening as being more pronounced at longer wavelengths (>0.9 μm) and larger phase angles (), such that a priori we should not expect a very pronounced effect (if any) in the visible-band PTF data.

Similar to the complication associated with comparing spin amplitudes from multiple oppositions (Section 3.1.1), an asteroid's mean color can potentially change if the spin axis varies with respect to our line of sight from year to year. Hence, we choose not to compare R-band and g-band phase-function fits from different oppositions. Aside from this constraint, we adopt the same two reliability selection criteria stated in Section 6.3, with a slight modification of requirement #2: here we allow four or more phase-angle bins of width , as opposed to the previous sections' five-bin requirement, because of the small sample size.

There are 92 asteroids with both R-band and g-band phase-function fits acquired during the same opposition that meet the above criteria. For each asteroid we difference the R-band G₁₂ value from the g-band G₁₂ value. The mean of this difference is , indicating (for the whole sample) no significant non-zero difference between the two bands' G₁₂ values. Likewise, for β, we compute a difference of , also consistent with zero difference between the bands.

Since these fits provide absolute magnitudes in each band (i.e., H_g and H_R), we compute the color for the 92-asteroid sample. Figure 20 shows that the distribution of this color is bimodal, suggesting that it is a viable proxy for taxonomy. This is further supported by the strong correlation between and the R-band G₁₂ value. No correlation is seen, however, between and the difference between the two bands' G₁₂ value or β values.

Zoom In Zoom Out Reset image size

The union of the color-index data (see Appendix) and photometric-index data (Section 6.3.1) provides significantly better taxonomic coverage of the PTF light curves (Figure 21). With this composite taxonomic information in hand, we can repeat the spin-amplitude analyses of Section 6.2 (Figure 17), this time considering the C-type and S-type groups separately. We define objects with one or both of the indices less than 0.25 as C type and greater than 0.75 as S type. We detail the resulting 1795-object taxonomically classified sample in Figure 22. There were 20 asteroids with conflicting color-based and photometric-based classifications that are not included in this 1795-object sample.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

The one-dimensional histogram in Figure 22 indicates that S-type asteroids dominates the smallest objects with data in PTF while C type dominates the largest. This reflects the fact that the survey's upper and lower sensitivity limits are defined in terms of absolute magnitude H (affected by albedo) rather than physical diameter, i.e., S-type asteroids larger than ∼50 km will tend to saturate the PTF detector, while C-type asteroids fainter than ∼5 km will usually fall below the detection limit. Adding to this effect is the fact that S types mostly occupy the inner main belt, where they are brighter by virtue of smaller heliocentric and geocentric distances, as compared to the usually more distant C types. While the two classes have similar representation in the sample (882 S types versus 913 C types), their true population ratio also affects the relative numbers.

The right-hand side plots in Figure 22 show rotation rate and amplitude versus diameter separately for the two taxonomic groups. Rather than plot a two-dimensional histogram as was done in Figure 17, for readability we here just plot the geometric mean and percentiles, computed by exactly the same running-bin method described in Section 6.2. The most prominent trend is that among asteroids, C types appear to rotate slower than S types and have larger amplitudes than S types. Assuming that both asteroid groups share the same mean angular momentum, the discrepancy could reflect the C types' ability to more efficiently redistribute material away from their spin axis, thereby increasing their moment of inertia (amplitude) while decreasing their angular rotation rate (i.e., a simple manifestation of conservation of angular momentum).

The above-stated assumption of a common mean angular momentum between C and S types is merely a simple case and is neither unique nor rigorously motivated. More careful consideration of, e.g., plausible ranges of internal tensile strengths of the two types could easily lead to more diverse scenarios wherein the two groups actually have different angular momenta and the observed spin-amplitude trends. As noted earlier (Section 6.2), large asteroids in general appear to have retained a significant primordial component in their spin distribution (Steinberg & Sari 2015); it is therefore important that differences in the origin of C types and S types (accretionary, temporal and/or spatial) be taken into account along with differences in collisional evolution and differing contributions from radiative forces like YORP. Simulations of the main belt's origin, such as the Grand Tack family of models (Walsh et al. 2011), should ultimately be modified to track particle spin evolution as well as orbits.

We also reproduce the period-versus-amplitude plot first shown in Figure 17, this time plotting separately the two taxonomic groups. The S types show a clearer cutoff at the 2 g cm⁻³ contour line, suggesting that they may in general be of greater bulk density than the C types, which show a softer boundary in this period-versus-amplitude space, the precise location of which appears to be somewhere between 1 and 2 g cm⁻³. Note that comparison to these density contours is only valid if the asteroids in consideration are held together mostly by self-gravity and approximated as fluids (as opposed to having significant internal cohesive or frictional resistance). These results are in general agreement with existing asteroid density estimates (Carry 2012 and references therein). Results from an independent analysis of a smaller, more densely sampled set of PTF asteroid light curves (Chang et al. in review; a study that follows closely the approach of Chang et al. 2014a) agree with the C type versus S type rotation rate discrepancy discussed here.

In order to compare the MPC-listed () magnitudes, which are in V band, with PTF's absolute magnitudes (, corresponding to the G-model fit), which are in either R or g band, we must first compute an approximate transformation from V band to each PTF band. While some transformations are given by Ofek et al. (2012a), we here prefer to empirically estimate these using actual asteroid photometry from both PTF and the MPC, rather than generating them from the more general transformations of Ofek et al. (2012a).

Figure 23 plots for asteroids whose PTF-derived is in the range . By restricting the comparison to objects with fitted values close to 0.15, we in principle select magnitudes for which the MPC's assumption is actually valid (none of the asteroids in Figure 23 have MPC-listed G values other than the default 0.15). Furthermore, we only consider (in Figure 23) asteroids with PTF data in at least three phase-angle bins of and either a reliable period or fitted amplitude less than 0.1 mag.

Zoom In Zoom Out Reset image size

Comparing the and magnitudes for this specific subset of asteroids, we obtain approximate transformations and . The uncertainties of 0.10 and 0.16 mag plausibly include a combination of the photometric calibration uncertainties of the MPC data (coming from a variety of surveys/observers), variation in H magnitude of a given asteroid between different oppositions (the MPC fits combine data possibly acquired at different viewing geometries), and the range of used in selecting the asteroids in this sample. Consideration of a range of values is equivalent to considering a range of asteroid colors (see the color-versus-G correlation seen in Figure 18). Hence, the uncertainties in these transformations also encompass the variation that might otherwise be formally fit in a color term for the transformations. Such a color term for R to V would almost certainly be less significant than that of g to V, as the former transformation is already zero within uncertainties. The larger uncertainty in the g to V transformation is likely attributable to both the smaller sample size and the fact that the V bandcenter is further displaced from g than from R, such that color variation has a more pronounced effect.

Given the above-computed transformations, and the fact that 89% of our fitted light curves are in R band, we proceed using only R-band light-curve fits, which we compare directly against MPC magnitudes (or formally, after applying the transformation of zero). A detail of the color dependence of the R to V transformation appears in the right plot of Figure 23; the mean transformation differs slightly between S and C types but not at a level comparable to the uncertainty in either.

In Figure 24 we show the relative error in the MPC absolute magnitudes as compared to the PTF magnitudes, for all 1630 light curves with sufficient phase-angle coverage in PTF (with the five-bin phase-angle criterion). These errors should reflect not only any discrepancy due to the different phase function models (PTF's G₁₂ versus MPC's G) but also variation in absolute photometric calibrations (within the MPC data internally and/or between the MPC and PTF data sets). The 0.1 mag uncertainty in the R to V band transformation has a prominent contribution to the errors shown here (the mean and 84th percentile of the errors expected from the 0.1 mag transformation uncertainty alone are shown as yellow dashed lines and assume p_V = 0.07). The green line (computed mean) and upper red line (84th percentile) indicate that the errors are ∼1% greater than those expected from the transformation uncertainty alone, though this increases slightly for the largest ( km) objects. Note that many of these largest asteroids are more frequently observed by programs other than the major sky surveys; these smaller facilities tend to use smaller aperture telescopes and different absolute calibration standards, which would contribute to the error.

Zoom In Zoom Out Reset image size

Instead of comparing just the fitted H magnitudes, for every light curve with a reliable PTF-fitted phase function we also compare the rms residual of all PTF data in that light curve with respect to both our G₁₂-fit-predicted R magnitude and the MPC (G = 0.15) predicted V magnitude. Our fit includes more fitted parameters and obviously should result in smaller residuals; Figure 25 shows that we see a factor of ∼3 smaller residuals in particular using the PTF fit. Note that if the 0.1 mag R-to-V transformation uncertainty were the only significant contributor to the MPC residuals, then their peak would instead be at ∼0.07 mag rather than ∼0.25 mag. Ignored rotational modulation and inaccurate phase functions move the MPC residuals distribution to higher rms values.

Zoom In Zoom Out Reset image size

The lower rms residuals afforded by the PTF light-curve model permit a more sensitive search for low-level transient activity (e.g., collisional events, cometary brightening) in these asteroids. For example, Cikota et al. (2014) perform a search for active main-belt asteroids using photometric residuals of all MPC data taken with respect to the MPC-predicted apparent V magnitudes. We currently are pursuing a similar analysis using these PTF light curves, as a follow-up to the morphology-based search already completed with PTF (Waszczak et al. 2013). A hybrid approach, wherein morphological measurements are made on stacked images of asteroids that have reliable light-curve fits, could further reveal this kind of subtle activity.

The purpose of this Appendix section is to introduce a one-dimensional color metric, based on data from seven different colorimetric asteroid surveys, which quantifies an asteroid's first-order visible-band color-based taxonomy as a number between 0 (C-type endmember) and 1 (S-type endmember). Our primary motivation for doing this is to enable a uniform comparison of PTF light-curve-derived parameters with color spanning from the brightest/largest objects (–9 mag, or –80 km diameters) down to PTF's detection limit for main-belt asteroids ( mag, or –4 km). Figure 26 panel A shows that the fraction of PTF light curves with color information increases by a factor of ∼3 among large asteroids when all seven surveys are considered, whereas for smaller objects the Sloan Digital Sky Survey's (SDSS; York et al. 2000; Ivezić et al. 2002; Parker et al. 2008) moving-object catalog provides essentially all of the color information.

Zoom In Zoom Out Reset image size

The seven surveys we use are described in Table 3. All of these surveys contain at least two independent color measurements, and when plotting their data in these two-dimensional spaces (or 2D subspaces defined by properly chosen principal components or spectral slope parameters), the first-order C-type and S-type clusters are in all cases prominently seen (Figure 26, panel B). To each such 2D color distribution we apply a two-dimensional fuzzy c-means (FCM) clustering algorithm (Bezdec 1981; Chiu 1994). For each survey data set, FCM iteratively solves for a specified number of cluster centers (in our case, two) in N dimensions (in our case, one dimension) by minimizing an objective function that adaptively weights each datum according to the robustness of its membership in a given cluster. The FCM output includes computed cluster centers and, for each datum, the probability that the datum belongs to each cluster (this being related to the datum's distance from each cluster center).

Table 3.  Asteroid Colorimetry Data Sets Used in Computing This Work's C/S Color Metric

Survey Name	References	Data Description	# Asteroids
UBV colors	Bowell et al. (1978)	U, B, and V broadband photometry acquired mostly at	902
	Tedesco (1995)	Lowell Observatory in the 1970 s with photomultiplier tubes.
Eight-Color Asteroid	Zellner et al. (1998)	Photometry in eight custom filters measured with photomultipliers at Catalina and
Survey (ECAS)	Zellner et al. (2009)	Steward Observatories. We compute and use the principal component color index	480
	⋯	PC#. Excludes objects with PC#1 error >0.3 mag.
24-Color	Chapman & Gaffey (1979, p. 655)	Photometry in 24 interference filters measured with photomultipliers at Mauna Kea.	262
Asteroid Survey	Chapman et al. (1993)	We compute and use the mean spectral reflectance slope and first principal component.
Small Main-belt	Xu et al. (1995)	CCD spectroscopy (0.4–1.0 μm, ) conducted mostly at Kitt Peak.
Asteroid Spectroscopic	Xu et al. (1996)	We compute and use the mean spectral reflectance slope and first principal component.	305
Survey (SMASS)	⋯	⋯
Small Main-belt	Bus & Binzel (2002)	CCD spectroscopy (0.4–1.0 μm, ) conducted at Kitt Peak.
Asteroid Spectroscopic	Bus & Binzel (2003)	We compute and use the mean spectral reflectance slope and first principal component.	1313
Survey II (SMASS-2)	⋯	⋯
Small Solar System	⋯	CCD spectroscopy (0.5–9.0 μm, ) conducted at ESO (La Silla).
Objects Spectroscopic	Lazarro et al. (2004)	We compute and use the mean spectral reflectance slope and first principal component.	730
Survey (S3OS2)	⋯	⋯
Sloan Digital Sky	Ivezić et al. (2002)	g,r,i, and z broadband CCD photometry acquired by SDSS from 1998–2009.
Survey (SDSS)	Parker et al. (2008)	Includes data in the Moving Object Catalog v4, supplemented with post-2007 detections	30,518
griz colors	Ivezić et al. (2010)	from SDSS DR10. We use the first principal component a* defined in the references.
	⋯	Excludes objects with a* error >0.05 mag or () error >0.1 mag.

Note. These data sets are visualized in Figure 26.

Download table as:  ASCIITypeset image

Table 4.  Parameters Describing PTF Light Curves with a Reliable Period or Phase Function. Byte-by-byte Description of file: ptf_asteroid_lc_parameters.txt

Bytes	Format	Units	Label	Explanations
1–4	I4	⋯	⋯	Light-curve ID numbera
6–11	I6	⋯	⋯	Asteroid number (IAU designation)
13–14	I2	yr	⋯	Last two digits of opposition year
16	I1	⋯	⋯	Photometric band: 1 = Gunn-g, 2 = Mould-R
18–20	I3	⋯	⋯	Number of observations in the light curve
22–26	F5.2	mag	⋯	Median apparent magnitude
28–37	F10.5	day		Time (MJD) of first observation
39–48	F10.5	day		Time (MJD) of final observation
50–54	F5.2	deg		Minimum-observed phase angle
56–60	F5.2	deg		Maximum-observed phase angle
62–63	I2	⋯	⋯	Number of sampled phase-angle bins of 3-deg width
65–68	F4.2	⋯	p	Reliability score from machine classifier: 0 = bad, 1 = good
70	I1	⋯	⋯	Manually assigned reliability flag: 0 = bad, 1 = good
72	I1	⋯	⋯	Period reliability flag: 0 = bad, 1 = good (product of two previous columns)
74–79	F6.3	mag	H	Absolute magnitude from G12 fit
81–85	F5.3	mag	⋯	Uncertainty in absolute magnitude from G12 fit
87–91	F5.3	⋯	G12	Phase-function parameter G12
92–98	F6.3	⋯	⋯	Uncertainty in
100–105	F6.3	⋯	G	Phase-function parameter G
107–113	F7.4	mag/deg	β	Phase-function parameter β
115–119	F6.3	mag	C	Phase-function parameter C
121–124	F4.2	mag	⋯	Amplitude from G12 fit (max − min)
126–134	F9.4	hr	P	Period from G12 fit
136–144	F9.4	hr	⋯	Period uncertainty from G12 fi
146–152	F7.4	mag	A11	Fourier coefficient from G12 fit
154–160	F7.4	mag	A12	Fourier coefficient from G12 fit
162–168	F7.4	mag	A21	Fourier coefficient from G12 fit
170–176	F7.4	mag	A22	Fourier coefficient from G12 fit
178–181	F4.2	⋯	⋯	Ratio of the two peak heights in folded rotation curvec
183–186	F4.2	⋯		Reduced chi-squared of the fit
188–192	F5.3	mag	⋯	"Cosmic error" (see Section 4.1)
194–198	F5.3	mag	⋯	Root-mean-square residual of observations w.r.t the fit
200–206	F7.3	hr	⋯	Reference period (from http://sbn.psi.edu/pds/resource/lc)
208–213	F6.2	km	D	Diameter derived from thermal IR datad
215–218	F4.2	km	⋯	Uncertainty in diameter
220–224	F5.3	⋯		Bond albedoe
226–231	F6.4	⋯	⋯	Uncertainty in bond albedo
233–236	F4.2	⋯	⋯	Color-based taxonomic index: 0 = C-type, 1 = S-type
238–241	F4.2	⋯	⋯	Photometry-based taxonomic index: 0 = C-type, 1 = S-type

Notes. ^aID number labels individual observations in Table 4. ^bSet to −1 if larger than the interval tested in grid search. ^cSet to 0 if there is only one maximum in the folded light curve. ^dReferences for the IR diameters are given in the text (Appendix). ^eBond albedo only computed for objects with reliable G₁₂ and available diameter.

Only a portion of this table is shown here to demonstrate its form and content. Machine-readable and Virtual Observatory (VOT) versions of the full table are available.

Download table as:  Machine-readable (MRT)Virtual Observatory (VOT)Typeset image

In the color-distribution plots of Figure 26, panel B (the plots with black backgrounds arranged diagonally), each pixel/bin is colorized according to the average cluster-membership probability of asteroids in that pixel. Blue indicates high probability of membership in cluster 1, while orange represents high probability of membership in cluster 2.

Our color index provides a more quantitative label than that offered by popular letter-based taxonomic systems (e.g., Bus et al. 2002 and references therein). Several such letter-based nomenclatures were in fact defined on the basis of one or more of these seven surveys, oftentimes using a method similar to the clustering technique we use here. We identify our blue cluster with C-type asteroids and our orange cluster with S-type asteroids, though we make this association purely for connection/compatibility with the literature. This is because our computed clusters have their own unique identity/definition, formally distinct from that given in any other work. Our clusters' definitions are nonetheless completely specified/reproducible by the FCM algorithm we used to compute them.

In reducing the taxonomic classification to a single number defined by the two most prominent groups (C and S types), we lose the ability to distinguish secondary classes like V types, D types, and so on. If such a sub-group is separated from both of the two main clusters, its members will be assigned membership scores of close to 0.5. For example, in the SDSS a* versus complex, the clearly seen V-type "tail" protruding down from the S-type cluster appears mostly green in color, reflecting its intermediate classification. Likewise for the less clearly seen D types, which in the SDSS plot lie above the S types and to the right of the C types (again in a green-colored region). The orders of magnitude lower numbers of such secondary types make them mostly irrelevant for the purpose of this analysis.

We compute the numerical uncertainty (variance) of a given asteroid's cluster-membership score in a particular survey by performing many bootstrapped trials wherein we first randomly perturb all data points by random numbers drawn from Gaussian distributions whose width are the quoted measurement (i.e., photometric) uncertainties in each of the two dimensions, and then repeat the FCM analysis on the perturbed data. The variance in each object's reported cluster probability is then computed after a large number of bootstrap trials.

Some asteroids appear in only one of the seven surveys; for such objects the color index is simply its cluster-membership score in that particular survey. For asteroids appearing in multiple surveys, we take the variance-weighted average of the multiple membership scores (and compute that composite score's variance by summing the component variances in inverse quadrature, as usual).

The many off-diagonal plots in Figure 26, panel B compare the cluster-membership scores of all asteroids appearing in all possible survey intersections. The number of asteroids in each survey (and in the intersection of each survey pair) appears above each plot (). The survey-pair distributions are 2D histograms where higher densities of data points correspond to black pixels/bins and low density or lack of data points is white. Evidently all possible survey combinations contain at least some asteroids (several share hundreds), and in all cases the individual taxonomic indices (on the horizontal and vertical axes) correlate strongly, confirming the consistency of the cluster membership between surveys.

In Figure 27 we illustrate some useful applications of this color index by comparing it with various asteroid surface observations. One of these quantities (SDSS a* color) was used in computing the color index, so its correlation with the clustering index is expected and thus confirmed.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

In the leftmost plot of Figure 27, asteroid photometry from GALEX²² (NUV band), compiled by A. Waszczak et al. (2015, in preparation), is normalized by the nominal G = 0.15 phase-model (Section 3.2.2) predicted brightness at the time of the GALEX observations, and the resulting NUV–V color evidently correlates with the visible color index. This indicates that asteroid reflectance slopes in the visible persist into the UV.

Figure 27 also plots our color index against the -band geometric albedo derived from WISE²³ observations obtained during its fully cryogenic mission. We only include asteroids that were detected in both of the thermal bands ( and ) and that therefore have a reliable diameter estimate. Use of this diameter in Equation (15) then permits estimation of the albedo, where the -band albedos use the corresponding WISE photometry (H in Equation (15) being replaced with the appropriate -band absolute magnitude).

The rightmost plot in Figure 27 shows our color index's relationship to a near-infrared color from the ground-based 2MASS (Skrutskie et al. 2006). Serendipitous asteroid detections were extracted from 2MASS by Sykes (2000, 2010) and include fluxes in J band (1.25 μm), H band (1.65 μm—not to be confused with the absolute visible magnitude H, used elsewhere in this work), and K band (2.17 μm).

Figure 28 plots our color index against proper orbital elements retrieved from the Asteroids Dynamic Site (AstDyS; Knežević & Milani 2012), revealing the distinct colors of dynamical families and the overall transition from S to C types with increasing semimajor axis. These are similar to the plots of Parker et al. (2008), which is not surprising given that the majority of the asteroids' color indices are based on SDSS data alone. Of the 32,5023 asteroids with a defined color index, there are 30,508 with proper orbital elements which are represented in Figure 28.

Similar to how we combined several surveys' colorimetric data in the previous section, here we compile thermal-infrared-derived diameter estimates from four surveys. Our aim is again to provide the largest possible sample for comparison with PTF-derived light curve data. Just as SDSS is the main contributor of colorimetry overall but suffers from incompleteness for large/bright asteroids, analogously WISE provides the vast majority of IR-based diameter measurements but levels off at ∼80% completeness at the bright end (Figure 29). We thus supplement WISE with diameter data from the Infrared Astronomical Satellite (IRAS; Matson et al. 1986; Tedesco et al. 2002), the Mid-Course Space Experiment (MSX; Tedesco et al. 2002), and AKARI Usui et al. (2011, 2014) compares several of these different data sets in terms of coverage and accuracy. As we did when defining the color index, asteroids occurring in multiple IR surveys are assigned the variance-weighted average diameter.

Regarding the WISE data in particular, we again use only those diameters that resulted from a thermal fit constrained by fluxes in all four WISE bands during the cryogenic mission. Furthermore, we use the latest (revised) diameter estimates published by Masiero et al. (2014), which adopted an improved thermal modeling technique first discussed by Grav et al. (2012).

The online version of this article includes two electronic tables containing the derived light-curve parameters and the individual photometric observations in each light curve. Tables 4 and 5 describe the columns and formatting of these tables, which include data on all reliable-period light curves, as well as those having amplitudes less than 0.1 mag and sampling in five or more 3-deg-wide phase-angle bins (which have reliable G₁₂ fits). Using these tables, one can produce plots of the PTF light curves we have analyzed in this work.