THE INFLUENTIAL EFFECT OF BLENDING, BUMP, CHANGING PERIOD, AND ECLIPSING CEPHEIDS ON THE LEAVITT LAW

The most common method to estimate parameters in a linear regression model is OLS, based on the Gauss–Markov theorem. The model is (in matrix notation)

$$\begin{eqnarray}&&{\boldsymbol{y}}={\boldsymbol{X}}{\boldsymbol{\beta }}+{\boldsymbol{\varepsilon }},\end{eqnarray} \tag{ 1 }$$

where ${\boldsymbol{y}}$ is the dependent variable vector, ${\boldsymbol{X}}$ is the independent variables matrix, ${\boldsymbol{\beta }}$ is the parameters vector, and ${\boldsymbol{\varepsilon }}$ is the errors vector.

This theorem states that the OLS estimators are the best linear unbiased estimators if the expected value of the errors vector is zero ($E({\boldsymbol{\varepsilon }})={\bf{0}}$), the variance of the errors is constant ($V({\boldsymbol{\varepsilon }})={\sigma }^{2}{\boldsymbol{I}}$), and the errors are uncorrelated. It can be shown (Montgomery et al. 2012) that OLS estimators ($\hat{{\boldsymbol{\beta }}}$) are obtained from

$$\begin{eqnarray}&&\hat{{\boldsymbol{\beta }}}={({{\boldsymbol{X}}}^{\prime }{\boldsymbol{X}})}^{-1}{{\boldsymbol{X}}}^{\prime }{\boldsymbol{y}}.\end{eqnarray} \tag{ 2 }$$

Notice that ${\bf{0}}$ is the null vector, ${\boldsymbol{I}}$ is the identity matrix, and ${{\boldsymbol{A}}}^{\prime }$ indicates the transpose of matrix ${\boldsymbol{A}}$.

Once the model is fitted its adequacy must be checked. That is: (i) the error term ${\boldsymbol{\varepsilon }}$ has a constant variance ${\sigma }^{2}{\boldsymbol{I}}$. (ii) The errors are uncorrelated. These items are the assumptions of the Gauss–Markov theorem. (iii) The errors are normally distributed. This is very important for testing the hypothesis and estimating confidence intervals. (iv) All the observations have approximately the same weight. This refers to the fact that it is not desirable that model estimators depend more on a few observations than on the majority of them. This could happen if some observations, called influential points, have a disproportionate impact in the OLS estimations. (v) There are no specification problems in the model and there is no structural break, i.e., the parameters are stable. This refers to the correct functional form of the model.

After fitting the model, a residual analysis is made to verify its adequacy. Since the residuals are the differences between the original observed values and their fits, they also measure the variability of the dependent variable that is not explained by the model. Therefore, any violations of the assumptions can be detected by analyzing model residuals. Mainly the analysis on residuals is looking for evidence of: (i) errors that come from a distribution with heavier tails than normal. Large departures from normality in the error distribution mean that the F- and t-tests are no longer valid, and neither are the estimations of the confidence interval. (ii) Heteroskedasticity, meaning that the errors do not have constant variance. (iii) Influential points. (iv) Specification problems. (v) Autocorrelated errors. In this study it is not necessary to check whether errors are correlated because the data are not time-dependent, i.e., they come from a cross-section in the statistical sense.

A normal probability quantil–quantil (Q–Q) plot is used for checking the error distribution. If all the residuals lie along a line, it means that they come from a normal distribution. For more details see Montgomery et al. (2012). This plot shows whether the distribution has heavier or lighter tails or whether it is skewed compared to the normal distribution.

To check the assumption that errors have constant variance (homoskedasticity) White's test (White 1980) is used in the following equation:

$$\begin{eqnarray}&&{\boldsymbol{y}}={\boldsymbol{X}}\hat{{\boldsymbol{\beta }}}+{\boldsymbol{e}},\end{eqnarray} \tag{ 3 }$$

where ${\boldsymbol{e}}$ is the residual from the linear regression defined as

$$\begin{eqnarray}&&\hat{{\boldsymbol{\varepsilon }}}={\boldsymbol{e}}={\boldsymbol{y}}-\hat{{\boldsymbol{y}}}.\end{eqnarray} \tag{ 4 }$$

This test does not make the assumption that errors come from a normal distribution. Its null hypothesis is that ${\sigma }_{i}^{2}={\sigma }^{2}$ for all i. Since the PL relation can be modeled as a simple linear regression, the test uses the auxiliary model given by

$$\begin{eqnarray}&&{e}_{i}^{2}={\alpha }_{0}+{\alpha }_{1}{x}_{i}+{\alpha }_{2}{x}_{i}^{2}+{u}_{i}\quad i=1,2,\ldots ,n,\end{eqnarray} \tag{ 5 }$$

where u_i is the error term and n is the number of stars. It can be shown that White's statistic nR² is asymptotically distributed as ${\chi }^{{}^{2}}$ with g – 1 degrees of freedom, where n is the number of observations and R² is the regression coefficient of determination of Equation (5). Because the auxiliary model (Equation (5)) has three parameters (g = 3), ${\chi }^{{}^{2}}$ has two degrees of freedom.

Outliers are points that have an unusual behavior. Let us imagine a scatter plot: outliers are data points that unusually are located out of the pattern or far from the data cloud. On the other hand, influential points are those data that have a disproportionate impact on OLS estimators (slope and zero point). As a result, the model estimators depend more on them than on the majority of data. Let us imagine the scatter plot again: a regression line is fitted to all points. Then we delete a point and fit the regression line once more. If the regression line obtained after removing the point changes a lot, then the deleted observation is an influential point. What is a lot? Statistics helps us in this decision, as we shall describe later. Now, let us imagine the scatter diagram with a point that is distant from the cloud of data points but located along the regression line. If we delete it, the regression line will not change. Therefore, the point is an outlier but not an influential point, because it affects the R² coefficient and the OLS estimator standard errors, but not the OLS regression estimators. In contrast, if the outlier affects the estimators of the regression coefficient it is also an influential point. Finally, let us imagine that the point is not an outlier; however, after deleting it, the new fitted regression line is very different from that obtained with it. This means that it is an influential point. In general, influential points may or may not be outliers.

Detection of influential points: in order to look for influential points Cook's distance, DFFITS, DFBETAS, and COVRATIO statistics are commonly used (Montgomery et al. 2012). The first of them is defined as follows:

$$\begin{eqnarray}&&{D}_{i}=\displaystyle \frac{{(\hat{{\boldsymbol{\beta }}}-{\hat{{\boldsymbol{\beta }}}}_{(i)})}^{\prime }({{\boldsymbol{X}}}^{\prime }{\boldsymbol{X}})(\hat{{\boldsymbol{\beta }}}-{\hat{{\boldsymbol{\beta }}}}_{(i)})}{p{\widehat{\sigma }}^{2}},\quad i=1,2,\ldots ,n.\end{eqnarray} \tag{ 6 }$$

Basically, it measures the square distance between the OLS estimate ${\boldsymbol{\beta }}$ based on all observations and the OLS estimate ${{\boldsymbol{\beta }}}_{(i)}$ obtained after deleting the ith observation, where p in Equation (6) is the number of parameters in Equation (3), i.e., the number of independent variables plus the zero point. Moreover, ${\widehat{\sigma }}^{2}$ is an estimate of the mean square error, defined as follows:

$$\begin{eqnarray}&&{\widehat{\sigma }}^{2}=\displaystyle \frac{{{\boldsymbol{e}}}^{\prime }{\boldsymbol{e}}}{n-p}.\end{eqnarray} \tag{ 7 }$$

Observations with large Cook's distance affect OLS estimates of ${\boldsymbol{\beta }}$. To know what is a large Cook's distance, D_i is compared to the ${F}_{0.5(p,n-p)}$ distribution (Montgomery et al. 2012). If the ith observation has ${D}_{i}\gt {F}_{0.5(p,n-p)}$, it is considered an influential point.

DFFITS measures how large the influence is of the ith observation on its own estimated ${\hat{y}}_{i}$. In other words, it shows how many standard deviations ${\hat{y}}_{i}$ changes due to the deletion of the ith observation. It is defined as follows:

$$\begin{eqnarray}&&{DFFITS}=\displaystyle \frac{{\hat{y}}_{i}-{\widehat{y}}_{(i)}}{{\hat{\sigma }}_{(i)}\sqrt{{h}_{{ii}}}},\quad i=1,2,\ldots ,n.\end{eqnarray} \tag{ 8 }$$

${\hat{y}}_{(i)}$ and ${\hat{\sigma }}_{(i)}$ are the fitted values and the standard deviation respectively, obtained without the ith observation. The term h_ii is the ith diagonal element of the matrix ${\boldsymbol{H}}$ that is defined as follows:

$$\begin{eqnarray}&&{\boldsymbol{H}}={\boldsymbol{X}}{({{\boldsymbol{X}}}^{\prime }{\boldsymbol{X}})}^{-1}{{\boldsymbol{X}}}^{\prime }.\end{eqnarray} \tag{ 9 }$$

The term h_ij is the amount of leverage that the ith observation has on the jth fitted value. A point is considered influential when $| {DFFITS}| \gt 2\sqrt{\tfrac{p}{n}}$. For complementary explanations and references about DFFITS and the ${\boldsymbol{H}}$ matrix see Montgomery et al. (2012) and Drapper & Smith (1998).

Another statistic that detects influential points is DFBETAS. It indicates how much effect the ith observation has on each ${\beta }_{j}$, measured in units of the standard deviation. In order to apply it, the coefficient estimators of ${\beta }_{j}$ obtained using all observations are compared with the coefficient estimators computed excluding the ith observation (${\beta }_{j(i)}$). As a result, a measure for each ${\beta }_{j}$ is obtained. The DFBETAS statistics is defined as follows:

$$\begin{eqnarray}&&{{DFBETAS}}_{j,(i)}=\displaystyle \frac{{\hat{\beta }}_{j}-{\hat{\beta }}_{j(i)}}{{\hat{\sigma }}_{(i)}\sqrt{{C}_{{jj}}}}\quad j=0,1,\ldots ,p\end{eqnarray} \tag{ 10 }$$

where ${\boldsymbol{C}}={({{\boldsymbol{X}}}^{\prime }{\boldsymbol{X}})}^{-1}$, so that C_jj is the jth diagonal element of the matrix ${\boldsymbol{C}}$ and ${\hat{\sigma }}_{(i)}$ is the square root of the regression mean square error fitted without the ith observation. An observation is considered an influential point if $| {{DFBETAS}}_{j,(i)}| \gt \tfrac{2}{\sqrt{n}}$.

The last statistic used to detect influential points is COVRATIO, which measures how the covariance matrix is affected by the ith observation. It compares the generalized variance of the parameter estimators obtained without the ith observation with that obtained using all observations.

The variance–covariance matrix of parameter estimators with all observations is

$$\begin{eqnarray}&&{Var}{(\hat{{\boldsymbol{\beta }}})={\hat{\sigma }}^{2}({{\boldsymbol{X}}}^{\prime }{\boldsymbol{X}})}^{-1}.\end{eqnarray} \tag{ 11 }$$

The generalized variance is the determinant of Equation (11):

$$\begin{eqnarray}&&{GV}(\hat{{\boldsymbol{\beta }}})=| {\hat{\sigma }}^{2}{({{\boldsymbol{X}}}^{\prime }{\boldsymbol{X}})}^{-1}| .\end{eqnarray} \tag{ 12 }$$

Given this, COVRATIO is defined as follows:

$$\begin{eqnarray}&&{{COVRATIO}}_{i}=\displaystyle \frac{| {\hat{\sigma }}_{(i)}^{2}{({{\boldsymbol{X}}}_{(i)}^{\prime }{{\boldsymbol{X}}}_{(i)})}^{-1}| }{| {\hat{\sigma }}^{2}{({{\boldsymbol{X}}}^{\prime }{\boldsymbol{X}})}^{-1}| }\quad i=1,2,\ldots ,n.\end{eqnarray} \tag{ 13 }$$

The ith observation is considered influential if its COVRATIO lies outside the interval $1-\tfrac{3p}{n}\lt {{COVRATIO}}_{i}\lt 1+\tfrac{3p}{n}$, where p is the number of parameters in Equation (3) and n is the total number of observations.

The LL reported in this work is obtained by rejecting from the sample the influential points and Cepheids exhibiting communalities, as explained later. However, two questions arise that require detailed answers: why is it important to identify these influential points and give them a statistical treatment? Why is it necessary to apply robust techniques to make a linear regression instead of applying a well-known algorithm such as sigma-clipping? In the next paragraphs, we answer these questions, illustrating the associated statistical problem.

Outliers are observations with anomalous behavior. When they are caused by human errors or instrumental failures, they can be recognized and excluded from the analysis. Moreover, they can be rejected from the analysis only if there are strong non-statistical reasons that support such a decision. In this work it is possible to reject them for astronomical reasons. However, the usual practice of deleting outliers, without making a further analysis to look for reasons that could explain their behavior, could lead to serious consequences in the precision of estimators, because they are adjusted artificially (Montgomery et al. 2012).

As a result, when there are no astronomical reasons to exclude those points and/or there are influential points, robust estimation methods are needed.

To establish that there is no evidence of specification errors in the model, we use tests to see whether a nonlinear combination of the fitted values is significant to explain the response variable. One of them is Ramsey's test, given by Equations (14)–(16) (Ramsey 1974). The other one is a variant proposed by Godfrey and Orme that uses only the model given by Equation (14) (Godfrey & Orme 1994). These auxiliary models and their respective hypotheses used in this work are the following:

$$\begin{eqnarray}&&{\boldsymbol{y}}={\boldsymbol{X}}\hat{{\boldsymbol{\beta }}}+{\alpha }_{1}{\hat{{\boldsymbol{y}}}}^{2}+{\boldsymbol{\varepsilon }}.\end{eqnarray} \tag{ 14 }$$

The null and alternative hypotheses are ${H}_{0}:{\alpha }_{1}=0$ and ${H}_{1}:{\alpha }_{1}\ne 0$

$$\begin{eqnarray}&&{\boldsymbol{y}}={\boldsymbol{X}}\hat{{\boldsymbol{\beta }}}+{\alpha }_{1}{\hat{{\boldsymbol{y}}}}^{2}+{\alpha }_{2}{\hat{{\boldsymbol{y}}}}^{3}+{\boldsymbol{\varepsilon }}.\end{eqnarray} \tag{ 15 }$$

The null and alternative hypotheses are ${H}_{0}:{\alpha }_{1}={\alpha }_{2}=0$ and ${H}_{1}:{\rm{a}}{\rm{t}}\;{\rm{least}}\;{\rm{one}}{\alpha }_{i}\ne 0;\quad i=1,2$

$$\begin{eqnarray}&&{\boldsymbol{y}}={\boldsymbol{X}}\hat{{\boldsymbol{\beta }}}+{\alpha }_{1}{\hat{{\boldsymbol{y}}}}^{2}+{\alpha }_{2}{\hat{{\boldsymbol{y}}}}^{3}+{\alpha }_{3}{\hat{{\boldsymbol{y}}}}^{4}+{\boldsymbol{\varepsilon }}.\end{eqnarray} \tag{ 16 }$$

The null and alternative hypotheses are ${H}_{0}:{\alpha }_{1}={\alpha }_{2}={\alpha }_{3}=0$ and H₁: at least one ${\alpha }_{i}\ne 0;$ $\quad i=1,2,3$.

If any of these tests are significant because there is evidence to reject the null hypothesis H₀, then there is enough evidence of misspecification. In fact, this means that a nonlinear combination of the fitted values is significant enough to explain the variability in the dependent variable, i.e., the functional form is incorrect.

Dummy variables are used to verify whether there are structural breaks. Before explaining what a dummy variable is, it is necessary to give a brief explanation about qualitative variables. Qualitative or categorical variables are those that do not have a natural scale of measurement. Two categorical variables are used in this study. The first variable indicates which data set the observation belongs to when data are divided into two sets, to make structural tests of the models. This variable has two categories: data set 1 and data set 2. The second variable is the name of the galaxy where the Cepheid is located, which also has two categories: LMC and SMC.

Dummy variables are used to indicate which category of a qualitative variable the observations belong to. Each dummy variable takes two possible values: 1 if the observation belongs to the category that the dummy variable represents, and 0 if not. If a qualitative variable has K categories, it can be represented by K dummy variables, one for each category.

However, only $K-1$ dummy variables are used in the regression model. Besides, it does not matter which $K-1$ dummy variables are used since the ${\boldsymbol{H}}$ matrix is the same independently of which dummy variable is omitted to fit the model. In fact, dummy variables representing the categories of a qualitative variable are linked to each other because when one of them takes a value of 1 for an observation, the others take 0 for the same observation. Therefore, when all the $K-1$ dummy variables in the model take a value of 0, they are representing the dummy variable that is not present in the model. For more details see Montgomery et al. (2012).

The process of verifying the existence of structural breaks is as follows: the data are divided into two data sets by dummy variables. Then, a model that includes the dummy variables is estimated to verify whether the estimated parameters associated with dummies are significant or not. If they are significant, the parameters are unstable, and thus there are structural breaks.

The P-value is used to contrast hypotheses. This value is the lowest level of significance that leads to rejection of the null hypothesis when it is true. Therefore, if the P-value is greater than the significance level, there is not enough evidence to reject the null hypothesis. In contrast, if the P-value is less than or equal to the significance level there is enough evidence to reject the null hypothesis and accept the alternative one, which means the test is significant.

Sometimes OLS assumptions are not satisfied; for example, errors have a distribution with heavier tails than the normal distribution and/or there are influential points. An alternative way to fit a linear model in scenarios of this kind is a robust regression, in which the residuals can be defined as follows:

$$\begin{eqnarray}&&{e}_{i}={y}_{i}-{{\boldsymbol{x}}}_{i}^{\prime }{\boldsymbol{\beta }}\quad i=1,2,\ldots ,n.\end{eqnarray} \tag{ 17 }$$

It is convenient to scale the residuals by using median absolute deviation:

$$\begin{eqnarray}&&S=\displaystyle \frac{\mathrm{median}| {e}_{i}-\mathrm{median}({e}_{i})| }{0.6745}.\end{eqnarray} \tag{ 18 }$$

Therefore the scaled residuals $({u}_{i})$ can be expressed as follows:

$$\begin{eqnarray}&&{u}_{i}=\displaystyle \frac{{y}_{i}-{{\boldsymbol{x}}}_{i}^{\prime }{\boldsymbol{\beta }}}{S}\quad i=1,2,\ldots ,n.\end{eqnarray} \tag{ 19 }$$

If scaled residuals are not dependent on each other, and they all have the same distribution f(u), the maximum likelihood ${\boldsymbol{\beta }}$ estimators are those that maximize the likelihood function:

$$\begin{eqnarray}&&L({\boldsymbol{\beta }})=\displaystyle \prod _{i=1}^{n}f\left(\displaystyle \frac{{y}_{i}-{{\boldsymbol{x}}}_{i}^{\prime }{\boldsymbol{\beta }}}{S}\right).\end{eqnarray} \tag{ 20 }$$

Drapper & Smith (1998) show that this condition is equivalent to minimizing Equation (21). Therefore, a class of robust estimator that minimizes this equation is called an M-estimator and the regression based on it is called M-regression. M comes from maximum likelihood since the function $\rho ({u}_{i})$ is related to Equation (20) for a proper choice of the error distribution.

Besides, $\rho ({u}_{i})$ weights the scaled residuals u_i in the sum

$$\begin{eqnarray}&&\displaystyle \sum _{i=1}^{n}\;\rho \left\{\displaystyle \frac{{y}_{i}-{{\boldsymbol{x}}}_{i}^{\prime }{\boldsymbol{\beta }}}{S}\right\}=\displaystyle \sum _{i=1}^{n}\rho ({u}_{i}).\end{eqnarray} \tag{ 21 }$$

There is more than one suggestion for $\rho ({u}_{i})$ in the literature; the function used in this work is Tukey's bi-square, with parameter equal to 6 (Kafadar 1983). For further information about this function and other choices see Drapper & Smith (1998). To minimize Equation (21) the first partial derivative of $\rho ({u}_{i})$ must be calculated:

$$\begin{eqnarray}\displaystyle \frac{\partial \rho }{\partial {\beta }_{j}} & = & \psi \left\{\displaystyle \frac{{y}_{i}-{{\boldsymbol{x}}}_{i}^{\prime }{\boldsymbol{\beta }}}{S}\right\}=\psi ({u}_{i})\\ i & = & 1,2,\ldots ,n,\quad j=0,1,2,\ldots ,p.\end{eqnarray} \tag{ 22 }$$

As a result, the following equation is obtained:

$$\begin{eqnarray}&&\displaystyle \sum _{i=1}^{n}\;{x}_{{ij}}\psi \left\{\displaystyle \frac{{y}_{i}-{{\boldsymbol{x}}}_{i}^{\prime }{\boldsymbol{\beta }}}{S}\right\}=0,\quad j=0,1,2,\ldots ,p\end{eqnarray} \tag{ 23 }$$

where p is the number of parameters in the model. Multiplying Equation (23) by ${u}_{i}/{u}_{i}$ gives

$$\begin{eqnarray}&&\displaystyle \sum _{i=1}^{n}\;{x}_{{ij}}{w}_{i}\left(\displaystyle \frac{{y}_{i}-{{\boldsymbol{x}}}_{i}^{\prime }{\boldsymbol{\beta }}}{S}\right)=0,\quad j=0,1,2,\ldots ,p\end{eqnarray} \tag{ 24 }$$

where ${w}_{i}=\psi ({u}_{i})/{u}_{i}$ ∀ ${y}_{i}\ne {{\boldsymbol{x}}}_{i}^{\prime }{\boldsymbol{\beta }}$, otherwise it is 1.

The determination of the M-estimator is an iterative process, thus the M-estimator for ${\boldsymbol{\beta }}$ is

$$\begin{eqnarray}&&{\hat{{\boldsymbol{\beta }}}}_{q+1}={({{\boldsymbol{X}}}^{\prime }{{\boldsymbol{W}}}_{q}{\boldsymbol{X}})}^{-1}{{\boldsymbol{X}}}^{\prime }{{\boldsymbol{W}}}_{q}{\boldsymbol{y}}.\end{eqnarray} \tag{ 25 }$$

This iterative process stops when a convergence criterion is reached, i.e., when the estimators change by less than a preselected amount in the ($q+1$)th iteration. Notice that Equation (25) is the weighted least-squares estimator where ${{\boldsymbol{W}}}_{q}$ is a diagonal n × n matrix of weights whose elements are w_i. See Drapper & Smith (1998) and Montgomery et al. (2012) for a wider discussion.

In the context of M-regression it is important to check that there are no bad influential points, i.e., observations that are outliers in the ${\boldsymbol{X}}$ space and influential points at the same time. In this scenario the MM-estimator (Yohai 1987) should be used because it has a high BDP. In a finite sample the BDP is the smallest fraction of contaminated data to cause a divergence of the estimator from the value that it would take if the data were not contaminated (Montgomery et al. 2012).

The MM-estimator combines the asymptotic efficiency of the M-estimator with the high BDP estimators such as S-estimators and Least Trimmed Squares. MM-regression is obtained in three steps. First, a model for high-BDP parameter estimators is calculated to compute the residuals of the model. Second, based on residuals computed in the previous step, an M-estimate of scale with high BDP is calculated. Finally, the model parameters are estimated using M-estimators and the scale estimation computed in the previous step.

Neither normal errors nor homoskedasticity are assumptions in M-regression or MM-regression. However, these regressions assume that errors are uncorrelated, which is not an issue when data come from a cross-section, meaning that data are not time-dependent.

Robust regressions are based on asymptotic results. This means they could lead to wrong results when they are used with small data sets. Montgomery et al. (2012) say that a small to moderate sample size could be less than 50 points. In our case, the smallest data set has about 400 observations. Hence, we have enough data to apply these techniques without problems. In particular, OGLE-IV has around twice and three times as many observations as OGLE-II for the SMC and the LMC, respectively.

In order to verify whether there are structural breaks (explained in Section 3.4) in model Equations (26)–(28), we fit the following model:

$$\begin{eqnarray}&&{\boldsymbol{y}}={\beta }_{1}+{\beta }_{2}{\boldsymbol{\delta }}+{\beta }_{3}{\boldsymbol{L}}{\boldsymbol{P}}+{\beta }_{4}{\boldsymbol{\delta }}{\boldsymbol{L}}{\boldsymbol{P}}+{\boldsymbol{\varepsilon }}.\end{eqnarray} \tag{ 32 }$$

Each data set is split into two parts using a dummy variable represented by δ, where

$$\begin{eqnarray}\delta =\left\{\begin{array}{ll}1 & {\rm{if}}\ P\ \leqslant \ 3.60\;(2.57)\;{\rm{day}}\\ 0 & {\rm{if}}\ P\ \gt \ 3.60\;(2.57)\;{\rm{day}}.\end{array}\right.\end{eqnarray} \tag{ 33 }$$

The LMC (SMC) data are split around P = 3.60 (2.57) days, in order to ensure that each data set has enough observations to avoid problems due to asymptotic results in robust regressions. As a result, the data sets have 240 (219) and 247 (235) points.

The hypotheses to contrast are ${H}_{0}:{\beta }_{2}={\beta }_{4}=0$ and ${H}_{1}:$ some ${\beta }_{j}\ne 0,\quad j=2,4$. The parameters ${\beta }_{2}$ and ${\beta }_{4}$ are associated with δ in order to detect structural breaks in the zero point $({\beta }_{2})$ and/or the slope $({\beta }_{4})$ of the LL; if these parameters are not significant there is not enough evidence of structural breaks. In contrast, if one or both of them are significant there is enough evidence of structural breaks.

The obtained LMC P-values are 0.886 for the V-band, 0.775 for the I-band, and 0.572 for the W_I index. Whereby, no evidence of structural breaks is found, hence the parameters obtained using M-estimators are stable at 5% significance.

We use the M-regression to estimate the LMC LL through the following equations:

$$\begin{eqnarray}&&{V}_{0}=17.005-2.683{LP}+e\end{eqnarray} \tag{ 34 }$$

$$\begin{eqnarray}&&{I}_{0}=16.549-2.930{LP}+e\end{eqnarray} \tag{ 35 }$$

$$\begin{eqnarray}&&{W}_{I}=15.839-3.315{LP}+e.\end{eqnarray} \tag{ 36 }$$

Models explain about 62%, 75%, and 87% of the variability in the V-band, I-band and W_I, respectively. The confidence intervals at 95% for each M-estimator for the LMC zero point $({\beta }_{1})$ and slope $({\beta }_{2})$ are reported in Table 3. It is worth mentioning that the estimators for zero point and slope obtained by OLS are included between the confidence intervals obtained by M-regression.

Table 3.  OGLE-II Leavitt Law for the LMC

	$\hat{{\beta }_{1}}$	$\hat{{\beta }_{2}}$	ZP	η
V0	17.005 ± 0.064	−2.683 ± 0.108	17.066 ± 0.021	−2.775 ± 0.031
I0	16.549 ± 0.044	−2.930 ± 0.074	16.593 ± 0.014	−2.977 ± 0.021
WI	15.839 ± 0.029	−3.315 ± 0.049	15.868 ± 0.008	−3.300 ± 0.011

Note. The reported β estimators are those given by Equations (34)–(36) and were obtained by using the M-regression. Confidence intervals are reported at 95%. The last columns report the values of zero point and slope obtained by Udalski (2000) using OLS.

Download table as:  ASCIITypeset image

To summarize, the robust tests have shown that there is no evidence to doubt the adequacy of the LMC models; as a result, they are acceptable.

For the SMC model the P-value obtained by the test of structural breaks is less than 0.001 for all bands, implying enough evidence to reject the null hypothesis and accept H₁ at 5% significance, i.e., there are structural breaks.

Because of the problems already detected in model Equations (26)–(28) for the SMC, further model analysis is needed using Equation (32). Therefore the individual significance of the parameters is analyzed. To do that, the following hypotheses are contrasted: ${H}_{0}:{\beta }_{j}=0$ and ${H}_{1}:{\beta }_{j}\ne 0,\quad j=2,4$. The individual robust test for ${\beta }_{2}$ shows it is significant at 5% significance, due to P-values of 0.035 in the V-band, 0.015 in the I-band, and 0.003 in the W_I index. In contrast, the same test for ${\beta }_{4}$ is not significant at 5% significance since the associated P-values are 0.566 in the V-band, 0.381 in the I-band, and 0.090 in the W_I index. There is not enough evidence to reject H₀, i.e., ${\beta }_{4}=0$. This means that the model has a different zero point for each data set but the same slope. Taking into account the specification problems and structural breaks detected in models (26)–(28), we fit the following model for SMC based on the previous individual significance analysis using (32)

$$\begin{eqnarray}&&{\boldsymbol{y}}={\beta }_{1}+{\beta }_{2}{\boldsymbol{\delta }}+{\beta }_{3}{\boldsymbol{L}}{\boldsymbol{P}}+{\boldsymbol{\varepsilon }}.\end{eqnarray} \tag{ 37 }$$

Using OGLE-III Cepheids, Subramanian & Subramaniam (2015) determined the inclination of the SMC galaxy's disk to be 64°. The geometrical distribution of the Cepheids across this disk could explain the behavior of the LL in the SMC, observed in Equation (37): linear regressions with the same slope and different zero points for each data set.

Model (37) is fitted by M-regression and no evidence of bad influential points is found. The robust Ramsey test is applied again to check whether there are specification problems in model (37). The equations and hypotheses of the test are

$$\begin{eqnarray}&&{\boldsymbol{y}}={\beta }_{1}+{\beta }_{2}{\boldsymbol{\delta }}+{\beta }_{3}{\boldsymbol{L}}{\boldsymbol{P}}+{\alpha }_{1}{\hat{{\boldsymbol{y}}}}^{2}+{\boldsymbol{\varepsilon }}.\end{eqnarray} \tag{ 38 }$$

The null and alternative hypotheses are ${H}_{0}:{\alpha }_{1}=0$ and ${H}_{1}:{\alpha }_{1}\ne 0$

$$\begin{eqnarray}&&{\boldsymbol{y}}={\beta }_{1}+{\beta }_{2}{\boldsymbol{\delta }}+{\beta }_{3}{\boldsymbol{L}}{\boldsymbol{P}}+{\alpha }_{1}{\hat{{\boldsymbol{y}}}}^{2}+{\alpha }_{2}{\hat{{\boldsymbol{y}}}}^{3}+{\boldsymbol{\varepsilon }}.\end{eqnarray} \tag{ 39 }$$

The null and alternative hypotheses are ${H}_{0}:{\alpha }_{1}={\alpha }_{2}=0$ and ${H}_{1}:{\rm{a}}{\rm{t}}\;{\rm{least}}\;{\rm{one}}\;{\rm{}}{\alpha }_{i}\ne 0;\quad i=1,2$

$$\begin{eqnarray}&&{\boldsymbol{y}}={\beta }_{1}+{\beta }_{2}{\boldsymbol{\delta }}+{\beta }_{3}{\boldsymbol{L}}{\boldsymbol{P}}+{\alpha }_{1}{\hat{{\boldsymbol{y}}}}^{2}+{\alpha }_{2}{\hat{{\boldsymbol{y}}}}^{3}+{\alpha }_{3}{\hat{{\boldsymbol{y}}}}^{4}+{\boldsymbol{\varepsilon }}.\end{eqnarray} \tag{ 40 }$$

The null and alternative hypotheses are ${H}_{0}:{\alpha }_{1}={\alpha }_{2}={\alpha }_{3}=0$ and ${H}_{1}:{\rm{a}}{\rm{t}}\;{\rm{least}}\;{\rm{one}}\;{\rm{}}{\alpha }_{i}\ne 0;\quad i=1,2,3$.

The test does not find enough evidence of specification problems because of the P-values obtained at 1% significance. In particular, for model (38) the SMC P-values are 0.933 for the V-band, 0.817 for the I-band, and 0.669 for the W_I index. For model (39), the SMC P-values are 0.370 for the V-band, 0.224 for the I-band, and 0.040 for the W_I index. For Equation (40), the SMC P-values are 0.503 for the V-band, 0.346 for the I-band, and 0.091 for the W_I index.

M-regression estimations for the SMC LL are shown in the following equations:

$$\begin{eqnarray}&&{V}_{0}=17.444+0.135\delta -2.551{LP}+e\end{eqnarray} \tag{ 41 }$$

$$\begin{eqnarray}&&{I}_{0}=17.019+0.107\delta -2.869{LP}+e\end{eqnarray} \tag{ 42 }$$

$$\begin{eqnarray}&&{W}_{I}=16.353+0.063\delta -3.342{LP}+e.\end{eqnarray} \tag{ 43 }$$

These models explain about 75%, 82%, and 87% of the variability in V, I, and the W_I index, respectively. The confidence intervals at 95% for each M-estimator in Equations (41)–(43) are reported in Table 4. Besides, it is interesting that the confidence intervals for M-estimators include the values obtained by OLS, although there are slight differences.

Table 4.  OGLE-II Leavitt Law for the SMC

	$\hat{{\beta }_{1}}$	$\hat{{\beta }_{2}}$	$\hat{{\beta }_{3}}$	ZP
V0	17.444 ± 0.097	0.135 ± 0.069	−2.551 ± 0.147	17.635 ± 0.031
I0	17.019 ± 0.074	0.107 ± 0.053	−2.869 ± 0.110	17.149 ± 0.025
WI	16.353 ± 0.052	0.063 ± 0.037	−3.342 ± 0.079	16.381 ± 0.016

Note. The reported β estimators are those given by Equations (41)–(43) and were obtained by using the M-regression. Confidence intervals are reported at 95%. The last column reports the zero-point value obtained by Udalski (2000) using OLS.

Download table as:  ASCIITypeset image

To summarize, robust tests show specification problems and structural break points in the SMC models (26)–(28); therefore, analysis leads to model (37) as a better choice.

The universality hypothesis of the LL implies that the slope of the linear regression observed in the same filter in different galaxies has the same value, and shows a negligible dependence on the metallicity (García-Varela et al. 2013, and references therein). In order to test the universality of the LL, a linear model is proposed to fit simultaneously the data from the LMC and SMC. This model is given by

$$\begin{eqnarray}&&{\boldsymbol{y}}={\beta }_{1}+{\beta }_{2}{\boldsymbol{\delta }}+{\beta }_{3}{\boldsymbol{L}}{\boldsymbol{P}}+{\beta }_{4}{\boldsymbol{\delta }}{\boldsymbol{L}}{\boldsymbol{P}}+{\boldsymbol{\varepsilon }}\end{eqnarray} \tag{ 44 }$$

which is fitted in V, I, and the W_I index after deleting the points previously discussed. The dummy variable δ indicates whether the observation is from the LMC or the SMC: $\delta =1$ if a Cepheid is from the SMC and $\delta =0$ otherwise. A total of 454 Cepheids from the SMC and 487 from the LMC are used to fit this model.

The model (44) is fitted by OLS. White's test finds enough evidence to reject the null hypothesis of homoskedasticity, at 5% significance, since the P-value is 0.017 in V, less than 0.001 in I, and less than 0.0001 in the W_I index. As happened before, OLS regression is inadequate due to heteroskedasticity problems and the presence of influential points. Besides, the distribution of residuals has a heavier tail than the normal distribution. Therefore, model (44) is fitted by M-regression using the bi-square function. No evidence of bad influential points is found; thus, the process to check the adequacy of the model continues.

The P-values of the slope test are less than 0.0001 in V, I, and for the W_I index, so the models are significant at 5% significance. The hypotheses contrasted are ${H}_{0}:{\beta }_{2}={\beta }_{3}\;={\beta }_{4}=0$ and ${H}_{1}:{\rm{s}}{\rm{ome}}\;{\rm{}}{\beta }_{j}\ne 0,\quad j=2,3,4$.

To test whether both galaxies have the same zero point and slope the following hypotheses are contrasted: ${H}_{0}:{\beta }_{2}={\beta }_{4}=0$ and ${H}_{1}:$ some ${\beta }_{j}\ne 0,\quad j=2,4$.

The P-values obtained by the robust test are less than 0.0001 in V, I, and for the W_I index; therefore, there is enough evidence to reject the null hypothesis at 5% significance. Hence, the individual significance of each parameter must be analyzed. To do that, the following hypotheses are contrasted again: ${H}_{0}:{\beta }_{j}=0$ and ${H}_{1}:{\beta }_{j}\ne 0,\quad j=2,4$.

The robust test of individual significance for ${\beta }_{2}$ shows that it is significant at 5% significance, since its associated P-values are less than 0.0001 in V, I, and the W_I index. The same test for ${\beta }_{4}$ is not significant at 5% significance because its associated P-values are 0.456 and 0.145 in V and I, respectively. This means that the model has different zero points, as we expect because the galaxies are at different distances, but it has the same slope. As a result, the final equations fitted by the M-regression in the V- and I-bands, after deleting ${\beta }_{4}$, are

$$\begin{eqnarray}&&{V}_{0}=17.037+0.554\delta -2.732{LP}+e\end{eqnarray} \tag{ 45 }$$

$$\begin{eqnarray}&&{I}_{0}=16.590+0.540\delta -2.995{LP}+e.\end{eqnarray} \tag{ 46 }$$

Confidence intervals are reported in Table 5. Models (45) and (46) make an adjustment to the zero point when Cepheids are from the SMC but have the same slope for both galaxies at 5% significance. Moreover, they explain about 78% and 81% of the variability in V and I, respectively.

Table 5.  Universal Leavitt Law Estimators

	$\hat{{\beta }_{1}}$	$\hat{{\beta }_{2}}$	$\hat{{\beta }_{3}}$
V0	17.037 ± 0.047	0.554 ± 0.030	−2.732 ± 0.074
I0	16.590 ± 0.034	0.540 ± 0.022	−2.995 ± 0.054

Note. The reported β estimators are those given by Equations (45) and (46), which were obtained by using the M-regression. Confidence intervals are reported at 95%.

Download table as:  ASCIITypeset image

The robust test of individual significance for ${\beta }_{4}$ in the W_I index is significant at 5% significance because its associated P-value is 0.0014. This implies that each galaxy has its own regression line with different zero point and slope, and it explains about 87% of the W_I variability. The following equation gives the model fitted for the W_I index:

$$\begin{eqnarray}&&{W}_{I}=15.841+0.594\delta -3.318{LP}-0.129\delta {LP}+e.\end{eqnarray} \tag{ 47 }$$

The confidence intervals at 95% for the estimated parameters in Equation (47) are $15.801\leqslant {\beta }_{1}\leqslant 15.880$ for zero point, $0.548\leqslant {\beta }_{2}\leqslant 0.639$ for zero-point adjustment when Cepheids are from the SMC, $-3.384\leqslant {\beta }_{3}\leqslant -3.251$ for slope, and $-0.209\leqslant {\beta }_{4}\leqslant -0.05$ for slope adjustment when Cepheids are from the SMC.

As a general result of this work, robust tests have shown that the LL of the LMC and SMC are universal in the V- and I-bands. There is enough evidence of two parallel regression lines, with the same slope at 5% significance but different zero points, as shown by models (45) and (46). However, for the W_I index, it has been shown that there is no universal LL: the LMC and SMC have two completely different regression lines, as shown by Equation (47).

Based on the experience obtained by fitting the models (26) and (27), they are fitted for all VI mean magnitudes of fundamental-mode Cepheids observed by the OGLE-II project. We also fit these models using the OGLE-IV fundamental-mode Cepheids, reported by Soszyński et al. (2015). It is worth mentioning that the published OGLE-IV mean magnitudes are not corrected for extinction. All the following models are fitted with the MM-regression because of its high BDP of 25%. The Least Trimmed Squares and Tukey function are used.

As no evidence is found of skewness in the models for both galaxies, we proceed with the analysis. Beginning with the LMC, model slope tests are significant, as can be seen in the P-values obtained: less than 0.0001 in the V- and I-bands (OGLE-II, OGLE-IV). Then, the test proposed by Godfrey & Orme (1994) is performed, giving P-values of 0.302 (OGLE-II) and 0.033 (OGLE-IV) in the V-band and 0.186 (OGLE-II) and 0.0003 (OGLE-IV) in the I-band; therefore, there is evidence of specification problems in the OGLE-IV I-band model at 1% significance. This could be due to blending problems and also to the lack of extinction correction in OGLE-IV data. To perform stability tests, the data are split by dummy variables into two data sets around P = 3.85 days (OGLE-II) and P = 3.5877 days (OGLE-IV), so as to have about half of observations in each set. No evidence of structural breaks is found, because the P-values obtained, by testing ${\beta }_{2}={\beta }_{4}=0$ versus the hypothesis that at least one of them is not 0 in Equation (32), are 0.601 (OGLE-II) and 0.453 (OGLE-IV) in the V-band, and 0.810 (OGLE-II) and 0.267 (OGLE-IV) in the I-band. Models explain about 61% and 68% of variability in the V- and I-bands respectively (OGLE-II), and 54% and 65% of variability in the V- and I-bands respectively (OGLE-IV). The estimated parameters of the models and their confidence intervals are reported in Table 6. Results for the I-band are shown only for information purposes because of their specification problems.

Table 6.  Optical Leavitt Law for the LMC

	OGLE-II		OGLE-IV		Soszyński et al.
	$\hat{{\beta }_{1}}$	$\hat{{\beta }_{2}}$	$\hat{{\beta }_{1}}$	$\hat{{\beta }_{2}}$	ZP	Slope
V	17.029 ± 0.041	−2.732 ± 0.063	17.436 ± 0.025	−2.701 ± 0.038	17.438 ± 0.012	−2.690 ± 0.018
I	16.553 ± 0.029	−2.947 ± 0.043	16.823 ± 0.018	−2.922 ± 0.029	16.822 ± 0.009	−2.911 ± 0.014

Note. The reported β estimators are those given by Equations (26) and (27). They were obtained by using the MM-regression. Confidence intervals are reported at 95%. The last columns report the values of zero point and the slope obtained by Soszyński et al. (2015) using the sigma-clipping algorithm and OLS over OGLE-IV data. The OGLE-II and OGLE-IV MM-estimators obtained in the V-band and their respective confidence intervals include at $1.0\sigma$ the OLS estimators reported by Udalski (2000) and Soszyński et al. (2015).

Download table as:  ASCIITypeset image

The OGLE-IV LMC data were also split into two around P = 10 days by Equation (32). No evidence of a structural break is found in the slope, because the P-value obtained is 0.705. However, the zero-point adjustment is significant. As a result, there is evidence of two parallel regression lines.

The SMC slope tests for models (26) and (27) are significant because the P-values obtained are less than 0.0001 for the V- and I-bands (OGLE-II, OGLE-IV). However, these models show problems. The test of Godfrey and Orme finds enough evidence of specification problems because the P-values obtained are less than 0.0001 for the V- and I-bands (OGLE-II, OGLE-IV). To perform stability tests, data are split around 2.0039 days (OGLE-II) and 1.9429 days (OGLE-IV) so as to have about half of observations in each data set. Enough evidence of structural breaks is found: P-values are less than 0.0001 in both bands (OGLE-II) when testing ${\beta }_{2}={\beta }_{4}=0$ versus the hypothesis that at least one of them is not 0 in Equation (32). Besides, when individual significance tests are performed for ${\beta }_{2}$ and ${\beta }_{4}$, the P-values obtained are less than 0.0001 for ${\beta }_{2}$ in both bands, and 0.0322 and 0.0082 for ${\beta }_{4}$ in the V-band and I-band respectively. As a result, in both bands there are two regression lines at 5% significance (OGLE-II).

Unlike OGLE-II models for the SMC, there is not enough evidence to support the idea that OGLE-IV model parameters for the slope are unstable with P-values 0.784 and 0.597 in the V- and I-bands respectively, while the zero-point adjustment, ${\beta }_{2}$ (Equation (32)), is still significant with a P-value less than 0.0001 in both bands. As a result, stability tests show that there are two parallel regression lines in both bands. New models are fitted based on the individual significance analysis using model (32) for the V- and I-bands. However, the test of Godfrey and Orme shows that specification problems persist. Nevertheless, for information purposes only, we report the estimated values of the slope at 95% confidence intervals for the OGLE-IV data: −2.780 ± 0.055 for the V-band and −3.024 ± 0.046 for the I-band.

Summarizing, the LMC models fitted with OGLE-II and OGLE-IV data are adequate except for the I-band model (OGLE-IV). The adequacy of the OGLE-II models could be explained because only $\sim 36 \%$ of Cepheids in the LMC have intermediate and large dispersion. Therefore, as Cepheids with small dispersion are the largest fraction, the exclusion of those with intermediate and large dispersion has no significant impact on the estimated slope value, because confidence interval for slope ($\hat{{\beta }_{2}}$), shown in Table 6, includes the slope estimation ($\hat{{\beta }_{2}}$) shown in Table 3. Moreover, LMC models fitted with all OGLE-IV data have confidence intervals for slope that contain estimations obtained by using the models fitted with all OGLE-II data (see $\hat{{\beta }_{2}}$ values in Table 6).

On the other hand, there is enough evidence to show that SMC models are inadequate. Unlike the LMC, the largest fraction of OGLE-II SMC data $(\sim 64 \% )$ are for Cepheids with intermediate and large dispersion; as a result, when models are fitted with all data, those Cepheids form the majority of the sample. The presence of a large fraction of Cepheids with intermediate and large dispersion in the SMC galaxy can be explained by the geometry of this galaxy, as discussed in Section 4.1. It is interesting that SMC models improve when they are fitted with OGLE-IV data because, despite specification problems, there is not enough evidence of instability of the slope parameter in both bands.