{'fraction_non_alphanumeric': 0.050517163843501724, 'fraction_numerical': 0.013441263178933693, 'mean_word_length': 4.164107857050703, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 1, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 23, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Doubly Robust Criterion for Causal Inference

The semiparametric estimation approach, which includes inverse-probability-weighted and doubly robust estimation using propensity scores, is a standard tool in causal inference, and it is rapidly being extended in various directions. On the other hand, although model selection is indispensable in statistical analysis, an information criterion for selecting an appropriate regression structure has just started to be developed. In this paper, based on the original definition of Akaike information criterion (AIC; \citealt{Aka73}), we derive an AIC-type criterion for propensity score analysis. Here, we define a risk function based on the Kullback-Leibler divergence as the cornerstone of the information criterion and treat a general causal inference model that is not necessarily a linear one. The causal effects to be estimated are those in the general population, such as the average treatment effect on the treated or the average treatment effect on the untreated. In light of the fact that this field attaches importance to doubly robust estimation, which allows either the model of the assignment variable or the model of the outcome variable to be wrong, we make the information criterion itself doubly robust so that either one can be wrong and it will still be an asymptotically unbiased estimator of the risk function. In simulation studies, we compare the derived criterion with an existing criterion obtained from a formal argument and confirm that the former outperforms the latter. Specifically, we check that the divergence between the estimated structure from the derived criterion and the true structure is clearly small in all simulation settings and that the probability of selecting the true or nearly true model is clearly higher. Real data analyses confirm that the results of variable selection using the two criteria differ significantly.


Introduction
Let us consider a fundamental setting in causal inference. It assumes that there are as many potential outcome variables as there are treatments, but that only the outcome variable corresponding to the assigned treatment is observed. There is confounding between the outcome and assignment, and if this is not taken into account and the estimation is done naively with the marginal likelihood, the estimates will have a bias. The bias can be avoided if the relationship between the outcome and confounding variables is modeled correctly, but a semiparametric approach without difficult modeling such as inverse-probability-weighted estimation (Robins et al. 1994) or doubly robust estimation (Scharfstein et al. 1999, Bang andRobins 2005) is often taken.
For example, let y is not observed, z i (∈ R) be the confounding variable. Furthermore, let us suppose a logistic model whose marginal probability function for the outcome variable is f (y where its regression structure is given by log{p (h) /(1 − p (h) )} = hθ, and whose conditional probability function for the as- i (see, for example, Hernán and Robins 2020).
Here, e (h) (z i ) is the propensity score introduced by Rosenbaum and Rubin (1983), and θ is the parameter relating to the causal effect. In this model, if we obtain the estimator of θ by i ; θ) despite the correlation between y i , it will have an asymptotic bias. The inverse-probability-weighted estimation is a method that uses the propensity scores to reproduce pseudo-complete data and gives an estimator by maximizing In this statistical problem, it is reasonable to treat the model selection of the regression structure relating to the causal effect as constant or quadratic, rather than linear hθ. However, for this basic selection, there are no reasonable information criteria. To be more specific, Platt et al. (2013) pioneered the information criterion for such models, proposing the use of inverseprobability-weighting for the goodness-of-fit term. While it is quite appropriate, their criterion uses the number of parameters of the regression structure as a penalty term, which will lead to a considerable underestimation of the bias correction. Baba et al. (2017) asymptotically evaluated the bias and corrected the penalty term. However, what they proposed is a C p -type criterion (Mallows 1973) that can only handle linear models, such as E(y (h) i ) = hθ, and cannot deal with the problem described above. In addition, it can not deal with the average treatment effect on the treated (ATT) or the average treatment effect on the untreated (ATU), which is often dealt with in causal inference, i.e., the causal effect when the target population is not necessarily the whole. In this paper, based on the original definition of Akaike information criterion (AIC; Akaike 1973), we derive an AIC-type criterion that overcomes these problems. If estimation using propensity scores is not used, the criterion proposed by Rolling and Yang (2014) may be a reasonable alternative; however, it should be noted that the semiparametric approach is important and is what we have assumed here.
When the propensity score e (h) (z) is unknown, a simple method is to assume some parametric function e (h) (z; α) as an approximation for it, obtain an estimatorα from {(t (h) i , z i ) : i ∈ {1, 2, . . . , N }}, and use e (h) (z;α) instead of e (h) (z). On the other hand, if we suppose that the model may be misspecified, we can instead assume some parametric function p (h) (y (h) | z; β) for the conditional distribution of the outcome variable given the confounding variable. Then, we can obtain an estimatorβ from {(t (h) i , y (h) i , z i ) : i ∈ {1, 2, . . . , N }}, and use p (h) (y (h) | z;β) to provide a doubly robust estimation; this is a standard method. If the modeling of either e (h) (z; α) or p (h) (y (h) | z; β) is correct, then θ can be estimated consistently. Although the information criterion given in Baba et al. (2017) covers doubly robust estimation, it is inadequate in that it is only derived when both modelings are correct. In this paper, we develop an information criterion that is an asymptotically unbiased estimator of an appropriate risk function as long as either modeling is correct, i.e., the information criterion itself is doubly robust.
Standard extensions of AIC include TIC (Takeuchi 1976) and GIC (Konishi and Kitagawa 1996); however, we would like to mention that our contribution does not fall within them. First, the penalty term in the AIC-type criterion is a certain asymptotic bias, so we usually have to consider convergence in mean, but to avoid unnecessary difficulties, we close the discussion with only weak convergence by using calculations peculiar to propensity score analysis. In the derivation of the doubly robust criterion, we use the asymptotic distribution of the doubly robust estimator under a setting that allows for model misspecification. To the best of our knowledge, there has been no effective use of that type of distribution. Moreover, the asymptotic bias depends both on the true distributions of the assignment and outcome variables, which may not be estimable in the setting of doubly robust estimation, and this is a difficulty not present in TIC and GIC; however, we propose a criterion that resolves this difficulty. It should also be emphasized that while the penalty terms of the traditional TIC and GIC are almost the same as twice the number of parameters, the penalty terms of our criterion tend to be much larger than that. This means that naively using the formal AIC can lead to the selection of significantly inappropriate models.
The organization of this paper is as follows. In Section 2, we describe a model, assumptions and a general causal effect in a setting where the target population is not necessarily the whole population, and introduce the inverse-probability-weighted estimator and the doubly robust estimator. In Section 3, we give a risk function based on the Kullback-Leibler divergence, which is naturally defined when considering such estimations, and define the AIC-type criterion by following the conventional derivation of AIC. Then, we asymptotically evaluate the penalty term of the information criterion for inverse-probability-weighted estimation. In Section 4, we derive the information criterion for doubly robust estimation, keeping in mind that one of the models of the assignment variable and the outcome variable may be misspecified. In Sections 5 and 6, we compare the performance of the proposed criterion with the criterion devised by Platt et al. (2013) through numerical experiments and an analysis based on real data. To explore the possibility of extending the derived information criterion, Section 7 tries to generalize divergence and weight functions, in order to treat, for example, a loss function robust to outliers and a covariate balancing propensity score. Finally, Section 8 summarizes our conclusions.

Model and assumption
In this paper, we treat a fundamental causal inference model, is an assignment variable that becomes 1 when the h-th treatment is is the probability function of y (h) given x (h) , and θ (∈ R p ) is the parameter used there (h ∈ {1, 2, . . . , H}). Note that y on the left-hand side is an observed outcome variable. Also, x (h) may contain some of the confounding variables, but for simplicity, the others are assumed to be non-random.
In this model, the (H − 1) potential outcome variables y (h) 's with t (h) = 0 are regarded as missing. In addition, since E(y (h) ) = E(y (h) | t (h) = 1) in general, naively estimating θ from the observed values alone will result in a bias. Here, we suppose that the confounding variable z (∈ R s ) for y (h) and t (h) is observed such that this bias can be removed. Moreover, we assume a weakly ignorable treatment assignment condition, which is intended to allow for the removal of this bias (Imbens 2000). We also assume positivity, P(t (h) = 1) > 0. There are N samples following this model, and the variables in the i-th sample are denoted with the subscript i. Letting u i = (y i , z i ) be the variables, the samples are assumed to be independent, i.e.
which naturally implies that y i ⊥ ⊥ y j .
The parameter of our interest to be estimated is θ in a population where the k-th treatment group is d (k) times larger than it actually is, that is, it can be regarded as the parameter for estimating a generalized average treatment effect (Kallus and Santacatterina 2022). If H = 2 and (d (1) , d (2) ) = (1, 0) or (d (1) , d (2) ) = (0, 1), θ is a parameter of the average treatment effect on the treated or the untreated. If H = 2 and (d (1) , d (2) ) = (1, 1), θ is a parameter of the average treatment effect of the whole sample. We denote the true value of the parameter by θ * , and where 0 p is the p-dimensional zero vector. We may actually should write (∂/∂θ) log f θ * (y i )| θ=θ * , but for simplicity, this notation is used here and hereafter.

Semiparametric estimation
If the relationship between the potential outcome variable y (h) and the confounding variable z can be modeled correctly, then the causal effect can be estimated consistently from the maximum likelihood method under the ignorable treatment assignment condition; however, this modeling is generally difficult. In recent years, a semiparametric approach using the propensity score E(t (h) | z), which does not necessarily require the above correct modeling, has often been used.
Hereafter, a model for E(t (h) | z) is denoted by e (h) (z; α) with a further parameter α (∈ R q ), and its true value is denoted by α * . Below, we discuss two typical estimation methods that are used in this approach.
The first is inverse-probability-weighted estimation (Robins et al. 1994). In this method, the missing values are pseudo-recovered by multiplying the observed values by the inverse of the propensity score as a weight, and then the usual estimation is implemented. Specifically, considering that the causal effect depends on d (k) , we define a weighted loss function using the weights w (h) (z; α * ) ≡ H k=1 d (k) e (k) (z; α * )/e (h) (z; α * ) and determine the inverse-probabilityweighted estimatorθ IPW by solving Conditional on z i , under the weakly ignorable treatment assignment condition, t and has the expectation e (h) (z i ; α * ), and then the left-hand side of (2) converges to the left-hand side of (1). This indicates that the inverse-probability-weighted estimator is consistent, i.e.θ IPW p → θ * . Moreover, by substitutingθ IPW into (2) and expanding the result, we obtain The above assumes that the propensity score is known, in other words, α * is known, but in general it is often unknown, in which case some estimatorα is substituted for α * . Forα, for example, we only have to construct the likelihood function for t (h) as a multinomial distribution with probability e (h) (z; α) and use the maximum likelihood estimator.
Although y (h) is correlated with z, the inverse-probability-weighted estimation does not directly use the information in z to estimate the expectation of y (h) . The doubly robust estimation (Scharfstein et al. 1999;Bang and Robins 2005) improves the inverse-probability-weighted estimation by doing so. Denoting a model for the conditional distribution of y (h) given z by p(y (h) | z; β) with a further parameter β (∈ R r ), and denoting its true value by β * , we take the expectation of log f (y (h) | x (h) ; θ) with p(y (h) | z; β) and write the conditional expectation as g (h) (x (h) , z; θ, β). In practice, some consistent estimatorβ is usually substituted for β * .
Specifically, if the propensity scores are unknown, we add to the left-hand side of (2), substitute estimators for α and β, and then find the doubly robust estimatorθ DR by solving with respect to θ such that it equals 0 p . This estimator not only improves the inverse-probability-weighted estimator but also is semiparametric locally efficient (Robins and Rotnitzky 1995). In addition, if either the propensity score or the conditional expectation is specified correctly, it is consistent and hence said to be doubly robust.

Risk function for causal inference
Before defining the risk function to derive the information criterion for causal inference, we will explain QIC w as proposed by Platt et al. (2013). QIC w is a criterion for missing data such as potential outcome variables, where the first term is the negative of twice the loss function used in the inverse-probability-weighted estimation, and the second term is twice the number of parameters. In the setting of Section 2.2, when α * is unknown, it is defined as When there are no weights t (h) i w (h) (z i ;α) in QIC w , QIC w is an asymptotic unbiased estimator of the risk function based on the Kullback-Leibler divergence. On the other hand, when weights are present, the variance of the first term increases, and the second term should be made larger.
Similarly to the usual AIC-type criteria, let (ỹ i.e., a random vector that is independently and identically distributed according to the distribution of (y i , z i ). Furthermore, letθ be an appropriate consistent estimator of θ with square root of n convergence, where we assume that we can write n 1/2 (θ − θ * ) = ξ + o P (1) using a random vector ξ which is O P (1) and whose expectation is 0. Here, we can define a risk such that the first term of (6) is naturally used as its naive estimator, where w (h) (z) is the limit of w (h) (z;α). Although some may be willing to use its true value, H k=1 d (k) E(t (k) | z)/E(t (h) | z), instead of the limit, here we consider the risk more matching toα, which we are actually using.
As we will see later, in fact, using the true value cannot provide a practical information criterion for doubly robust estimation. Note that when α * is known, we use w (h) (z; α * ) as w (h) (z). This risk function can be regarded to be based on the Kullback-Leibler divergence between the true and estimated distributions in a population where the k-th treatment group is d (k) times larger than it actually is. Also, as can be seen from (2), the estimation aims to minimize this risk function, and in that sense, it is a natural quantity. Since the first term in (6) is an evaluation of this risk function using the same data as those for the estimator, it tends to be smaller than the true value. Therefore, denoting the bias as we define b limit as the weak limit of the quantity in the expectation in (8), i.e., the random variable to which the quantity converges in distribution, and use E(b limit ) as the asymptotic bias for the correction. If (8) is expanded with respect toθ around θ * , we get The expectation of the second term is divided into the expectation of ξ and the expectation of the other terms; the former is 0 p , so we obtain When this asymptotic bias contains unknown parameters, they are replaced by their consistent estimators, as is done in deriving the usual AIC-type information criterion, and the evaluated asymptotic bias is added to the first term of (6) to construct the information criterion.
This construction of the information criterion is based on Baba et al. (2017), but the same kind of construction is used in Wallace et al. (2019), so we would like to mention the differences between the two. Wallace et al. (2019) proposes a valuable information criterion that is also justified in causal inference, and actually the main focus is on its application to g-estimation. On the other hand, the derivation is intuitive, and the information criterion is defined by substituting an estimator into the expectation of quasi-log-likelihood. Specifically, letting Q(·) be the quasilog-likelihood,ψ be the estimator, and ψ † be the limit of the estimator, then −2Q(ψ) is used instead of E{−2Q(ψ † )} without considering the bias. In this paper, while being conscious of the difference between convergence in mean and weak convergence, we use calculations peculiar to propensity score analysis based on the ignorable treatment assignment condition, and construct the information criterion in a way that may seem somewhat complicated. As will become clear later, while the penalty term of the criterion in Wallace et al. (2019) for g-estimation should be close to that of AIC, the penalty term of our criterion for estimation based on inverse-probabilityweighting tends to be much larger than that of AIC.

Case of known propensity score
Since the error of the inverse-probability-weighted estimator can be written as (3), we hereafter evaluate The last equality is derived from (1). Therefore, the expectation in the right-hand side of (10) becomes the trace of Therefore, by defining we arrive at the following theorem.
Theorem 1. Suppose that A(θ, α) and B(θ, α) are defined as in (4) and (13). Then, the asymptotic bias of the information criterion for the inverse-probability-weighted estimation is given by when the propensity score is known.
On the basis of this result and the fact thatθ IPW is a consistent estimator of θ, we propose as an information criterion for the inverse-probability-weighted estimation when the propensity score is known, whereÂ(θ, α) andB(θ, α) are empirical version of A(θ, α) and B(θ, α), from the definition ofB(θ, α) is the same asÂ(θ, α). It means that if our target of estimation is close to the average treatment effect of the whole sample, the IPWIC penalty term will be much larger than twice the number of parameters.

Case of unknown propensity score
If the parameter α in the propensity score e (h) (z; α) is unknown, then we only have to and the error ofα is expressed aŝ From this and (3), the error of the inverse-probability-weighted estimator iŝ where Since the error can be written like this, on the basis of (9), we hereafter evaluate When i = j, this expectation is divided into one for i and one for j because of their independence, and the sum of the expectations of the terms with i taken over h is 0 p from (11). Therefore, we only have to evaluate the terms with i = j. When k = h, the components are 0 because On the other hand, when k = h, the first term is the same as (12). The next term By defining we arrive at the following theorem.
On the basis of this result and the fact that (θ IPW ,α) is a consistent estimator of (θ * , α * ), we propose as an information criterion for the inverse-probability-weighted estimation when the propensity score is unknown, whereÎ 1 (α),Λ 1 (θ, α) andΛ 2 (θ, α) are empirical version of I 1 (α), Λ 1 (θ, α) and Λ 2 (θ, α), respectively. Theorems 1 and 2 indicate that the penalty for an unknown propensity score tends to be smaller than that for a known propensity score. This is consistent with the fact that the asymptotic variance for the inverse-probability-weighted estimator becomes smaller if the propensity score is estimated even when it is known (see, e.g., Henmi and Eguchi 2004).

Doubly robust criterion
The doubly robust estimator is an estimator that is consistent even if either the model for the assignment variable conditional on the confounding variable or the model for the outcome variable conditional on the confounding variable is misspecified. In the narrow setting of a linear model and a basic average treatment effect, Baba et al. (2017) proposed a C p -type criterion for the doubly robust estimation, but derived it under the assumption that both models are correct, so it has no validity when only one of them is correct. In this section, we aim to derive a doubly robust criterion that is an asymptotically unbiased estimator of the risk function in (7). Note that in the setting of doubly robust estimation, the propensity score is unknown; i.e., α * is unknown.
Then, the doubly robust estimating equation is Let us denote the limit of the solution (θ DR ,α,β) of this estimating equation by (θ † , α † , β † ). If one of the models is correct, θ † = θ * (we will write θ * from now on), but not necessarily α † = α * or β † = β * . We perform a Taylor expansion on the estimating equation with the estimator substituted into it and develop a statistical asymptotic theory similar to the conventional one.
The derivative of the left-hand side of the estimating equation is asymptotically the sum of and if we divide it by −N , it converges in probability to from the law of large numbers. Here, I 1 (·) and I 2 (·) are the Fisher information matrices for α and β, respectively. and O is a zero matrix. The terms in the (1,1) and if either model is correct, it becomes A(θ * , α † ) on the basis of the definition of (4). This is because, if the model for the assignment variable is correct, then E{t (h) , whereas if the model for the outcome variable is correct, then and it can be seen that the inverse of (20) is Note that in (20), the terms in the (1,2) block are if the model for the outcome variable is correct, and the terms in the (1,3) block are if the model for the assignment variable is correct. This means that Using the representation of (23), we obtain the following lemma from conventional statistical asymptotic theory.
Lemma. Suppose that A(θ, α), m(u; θ, α, β), C 1 (θ, α, β), and C 2 (θ, α, β) are defined as in Here, at least one of α † or β † will be the true α * or β * . In particular, when both models are correct, the error is expressed aŝ As an expansion formula forθ DR , (25) is used to show the local asymptotic efficiency, but (24) is used here to derive an information criterion while preserving double robustness, which is unexampled to the best of our knowledge. In order to derive the information criterion using this result and (9), we will evaluate First, let us evaluate the sum of expectations in the first term with respect to h. When i = j, the expectation is divided into one for i and one for j. If either the model for the assignment variable or the model for the outcome variable is correct, the latter E{m(u j ; θ * , α † , β † )} is 0 q .
Therefore, we only have to consider the case of i = j; letting it can be seen that the sum becomes Next, let us evaluate the sum of the expectations in the second term with respect to h. When i = j, the expectation for j is H k=1 E{t (k) j ∂ log e (k) (z j ; α † )/∂α} = 0. Therefore, we only have to consider the case of i = j; it can be seen that the sum becomes Finally, let us evaluate the sum of the expectations in the third term with respect to h. When j | z j ; β † )}/∂β] = 0. Therefore, we only have to consider that case of i = j, and it can be seen that the sum becomes From the above, we arrive at the following theorem.
are defined as in (4), (13), (26), (27) and (28). Then, if either the model for the outcome variable conditional on the confounding variable or the model for the assignment variable conditional on the confounding variable is correct, the asymptotic bias of the information criterion for the doubly robust estimation is given by When the former model is correct, D 2 (θ * , α † , β † ) becomes 0; when the latter model is correct, Since this asymptotic bias depends on both the true distribution of the outcome variable conditional on the confounding variable and the true distribution of the assignment variable conditional on the confounding variable, one may think that Theorem 3 is not practical. However, instead, an empirical estimation can providê which are consistent for D 1 (θ * , α † , β † ), D 2 (θ * , α † , β † ), and D 3 (θ * , α † , β † ), respectively. Note (7), we would not have these consistent estimators, so the asymptotic bias in that case is not practical. From the above, we propose as an information criterion for the doubly robust estimation.

Estimation of various average treatment effects in continuous outcomes
Let us consider a model with a treatment group represented by h = 1 and a control group represented by h = 2, i.e. H = 2, where the outcome variable y (h) (∈ R) is observed for either of them. As the explanatory variables that compose the regression structure, we suppose that x = (x 1 , x 2 , x 3 , x 4 ) ′ is distributed according to N(0 4 , I 4 ), and x (1) = (x 1 , x 2 , 0, 0) ′ and x (2) = (0, 0, x 3 , x 4 ) ′ , which are different for each group. The assignment variable that becomes 1 when the sample is in the treatment group is represented by t (1) , and the assignment variable that becomes 1 when the sample is in the control group is represented by t (2) . We assume that y (h) and t (h) are correlated and that the confounding variable z = (z 1 , z 2 ) ′ which can explain the correlation follows a Gaussian distribution N(0 2 , I 2 ) independently of x.
As the model for t (h) , we will use a logit model in which z is the explanatory variable.
The target of the estimation is the treatment effect with (d (1) , d (2)  First, let us examine whether the IPWIC penalty term in (18)  reasonably well in all cases. As a matter of fact, when we consider the ATT in two groups, the effect of the inverse probability on the bias is reduced, so that even the penalty term of QIC w gives a reasonably good approximation. Actually, when d = 0, the Monte Carlo evaluation (MCE) is close to 4; here, it cannot be concluded whether the penalty term of QIC w is superior or inferior to the penalty term of the proposed criterion. On the other hand, if we consider the ATE, the inverse probability effect appears, so the MCE at d = 1 is larger than 4 and the proposed penalty term is superior. As the number of groups increases, the effect of the inverse probability becomes more pronounced, and the evaluation of the proposed penalty term becomes more meaningful even when the ATT is considered. Note that this subsection only treats the simple two-group case; the multi-group case is treated in the next subsection.
Tables 2 and 3 compare the results of model selection for the settings where the bias evaluation of the proposed criterion works well, with and without model misspecification. Since the propensity score is unknown, we consider IPWIC in (18) and DRIC in (29) to be the proposed criteria and QIC IPW w using inverse-probability-weighted estimation in QIC w and QIC DR w using doubly robust estimation in QIC w to be the comparison targets. As models, we examine the nine types described above and select the optimal one for each criterion. As the main index to measure the goodness of the criteria, we use the empirical estimator of the risk function in (7) minus a constant independent of the model selection, specifically, the average of 3000 Comparing the existing criteria with the corresponding proposed criteria, i.e., QIC DR w and DRIC, and QIC IPW w and IPWIC, we see that the proposed criteria are superior in all cases, including the case of d * = 0 where QIC w should be somewhat appropriate. According to the selection probabilities, we can see that when d * = 0.5 or d * = 1, the existing criteria overselect Model 2 or Model 4 and lead to overfitting, and that the effect of underestimating the penalty term is strongly apparent. Note that model misspecification does not degrade performance as the theory suggests and that DRIC is always superior to IPWIC.
Of course, it holds that P(y (h) = 1 | z; θ, β) = 1 − P(y (h) = 0 | z; θ, β). Letting ǫ be a random effect that follows N(0, 1 − β ′ β) independently of z, the model of y (h) is a logit model with the regression function as θ 1 + θ 2 h + θ 3 h 2 + θ 4 h 3 + β ′ β + ǫ. We suppose which marginalizes z, as the probability function of y (h) in the model. In generating data according to P(y (h) = 0 | z; θ, β) above, the true value of the parameter to be used is the one that satisfies θ * 1 + θ * 2 h + θ * 3 h 2 + θ * 4 h 3 = 0.5 + θ * (h − 1) in the case of a linear polynomial, θ * 1 + θ * 2 h + θ * 3 h 2 + θ * 4 h 3 = 0.5 + 0.2(h − 1) + θ * (h − 1) 2 in the case of a quadratic polynomial, and in both cases (α * 1 , α * 2 , α * 3 , α * 4 , α * 5 ) = (0.2, −0.15, −0.25, 0.4, α * ) and (β * 1 , β * 2 ) = (0.2, β * ). We will set 0.05, 0.01 or 0.005 for θ * , and 0 or 0.2 for α * and β * . As candidate models, let us consider a 0th-order polynomial (constant) model with θ 2 = θ 3 = θ 4 = α 5 = β 2 = 0, a 1st-order polynomial (linear) model with θ 3 = θ 4 = α 5 = β 2 = 0, a 2nd-order polynomial (quadratic) model with θ 4 = α 5 = β 2 = 0, and a 3rd-order polynomial (cubic) model with α 5 = β 2 = 0. That is, in all of these models, if α * = 0, then the model for the assignment variable has been misspecified, and if β * = 0, then the model for the outcome variable has been misspecified. The target of the estimation is the average treatment effect of the whole sample 1, 1, 1). The sample size is N = 100 or N = 200, and the number of repetitions is 3000. Now let us examine whether the IPWIC penalty term in (18) and the DRIC penalty term in (29) can approximate the true bias represented by (8) in Table 4. In this table, when the true model is linear, the assumed model is also linear, and when the true model is quadratic, the assumed model is also quadratic. That is, in QIC w , the penalty term is 2 × 2 = 4 when the model is linear, and 2 × 3 = 6 when the model is quadratic. The table confirms that the approximation works reasonably well in all cases. At least, they are much closer to the true bias than 4 or 6 in the evaluation of QIC w .  Table 5 compares the results of the model selection for the settings where the bias evaluation of the proposed criterion works well. Since the propensity score is unknown, we consider IPWIC in (18) and DRIC in (29) to be the proposed criteria and QIC IPW w using inverse-probabilityweighted estimation in QIC w and QIC DR w using doubly robust estimation in QIC w as the comparison targets. As models, we examine 0th-, 1st-, 2nd-and 3rd-order polynomials in h as described above and select the optimal one for each criterion. As the main index to measure the goodness of the criteria, the evaluated value of risk in (30) is used as in Section 5.1, and also the selection probability is checked as a reference index. Comparing the existing criteria with the corresponding proposed criteria, i.e., QIC DR w and DRIC, and QIC IPW w and IPWIC, we find that the proposed criteria are clearly superior in all cases. According to the selection probabilities, the existing criteria are overfitting by selecting too large an order, and the effect of underestimating the penalty term is rather apparent. Especially when θ * is small, one might be concerned that the proposed criterion chooses a first-order model that is not true with considerable probability. However, if the contribution of the quadratic term is small, selecting the true model is not necessarily optimal because of the principle of parsimony. Actually, we can confirm that the risk of the proposed criterion is small even if non-true models are selected to some extent.

Real data analysis
The LaLonde dataset is treated in LaLonde (1986) and is included in the R package Matching.
The group that took the U.S. job training program in 1976 is denoted as t = 1 and the group that did not take the program is denoted as t = 0. The difference in annual income in 1978 after the training for each group is estimated as the average treatment effect on the treated (ATT) or the average treatment effect of the whole sample (ATE). The confounding variables are age The MineThatData dataset is published by Hillstrom (2008) and contains 12 attributes, such as purchase amount and district classification code, of 64000 customers. The purpose here is to estimate the causal effect of an email sent to male customers on the probability of the men purchasing a product. Since the sample size of 64000 is too large a computational load, we use a sample of N = 25508, where the district code (zip code) is not urban, and divide the sample into two groups: a group to which the email was delivered to the men (t = 1) and a group to which it was not sent (t = 0). The confounding variables are the number of months since the last purchase (recency), the amount purchased in the last year (history), whether the user purchased men's products in the last year (men's), whether the user purchased women's products in the last year (women's), whether the user became a new user in the last 12 months (newbie), i.e. p = 5, and the outcome variable is whether the user made a purchase within two weeks of receiving the email (conversion).
For these datasets, we use the same model as in Sections 5.1 and 5.2, and assume that the regression structure for each latent variable can be written as a linear sum of confounding variables. The confounding variables are then selected by using the step-up procedure using each criterion; the regression coefficients of the selected variables are shown in Table 6. For the LaLonde dataset, we examined both the ATT-targeted DRIC and QIC w and the ATE-targeted DRIC and QIC w ; while the former showed some differences between the two criteria, the latter showed larger differences. For the MineThatData dataset, the proposed criterion selects three variables less than the existing criterion, from which we can confirm that the same trend as in the numerical experiment appears in the real data analysis. Since we do not know the true structure of the real data, it is impossible to judge the superiority or inferiority of the two criteria; however, we can see that the difference between them is quite large.

Generalization of divergence
While the doubly robust estimation is literally robust against misspecification of the model, it is not robust against outliers because it uses the log-likelihood as the loss function, and Table 6: Estimates for the real data given by each criterion. Variables with an estimated value of 0 imply that they were not selected in the model selection. A value of 0 attached to the criteria means that it is for ATT, and a value of 1 means that it is for ATE. The exponential notation is used in the figures, and "E" is omitted.
Then, it can be seen that all the content on doubly robust estimation in this paper applies after replacing log f with ζ and g (h) with η (h) . For example, we can get the doubly robust estimating equation by adding (5) to the left-hand side of (2). When the model for the propensity score is correct, the expectation of the quantity in the curly brackets in (5) is 0, so the estimation equation converges to (31). Also, when the model for the outcome variable is correct, the expectations of t (h) w (h) (z; α † )ζ(y (h) | x (h) ; θ * ) and t (h) w (h) (z; α † )η (h) (x (h) , z; θ * , β † ) are equal, so asymptotically, these terms cancel out and the estimating equation still converges to (31).
This means thatθ DR converges in probability to θ * .
The content in Section 4 is only valid when either α † or β † is the true value, so one may be concerned that θ * is not necessarily the true value, but there is no part which is affected by that. Accordingly, Theorem 3 holds. The true structure appears in the construction ofD 1 ,D 2 andD 3 , but again, whether θ * is the true value or not has no effect. Actually,D 1 ,D 2 andD 3 are constructible as long as E(t (h) | z) and E{∂ log f (y (h) | x (h) ; θ)/∂θ ′ | z} appear one by one in product form. As a result, we have devised a criterion which replace log f with ζ and g (h) with η (h) in (29).

Generalization of weight function
Up to this point, the weight function w (h) (z; α) for t (h) has been H k=1 d (k) e (k) (z; α)/e (h) (z; α) for both the derivation of the estimator and the derivation of the information criterion. In fact, it is possible to construct similar information criteria as in Sections 3 and 4 even if we change the weight function in each of the derivations; two examples of this will be presented below.
The first example focuses on the true structure for the observed data, as was addressed in Baba et al. (2017). Specifically, it corresponds to changing only the weight function used in deriving the information criterion to H k=1 d (k) e (k) (z; α). In this case, for example, for the inverse-probability-weighted estimation with a known propensity score, Theorem 1 holds if we remove 1/e (h) (z; α) from the definition of B(θ, α) in (13). Furthermore, when d (1) = d (2) = · · · = d (H) = 1, which is a setting that considers the average causal effect of the whole sample, it holds that A(θ, α) = B(θ, α) because H k=1 d (k) e (k) (z; α) = H k=1 e (k) (z; α) = 1. Then, the asymptotic bias in Theorem 1 is calculated to be −2p. In other words, although Platt et al. (2013) does not consider this weight function, the information criterion in this case can be written as (6).
The second example uses covariate balancing (Imai and Ratkovic 2014), in which the weight function used to derive the estimator is changed. The basic idea of covariate balancing is not to use the maximum likelihood estimator based on the modeling of e (h) (z; α) for providingα but rather to use the information in E{(∂/∂θ) log f (y (h) | x (h) ; θ) | z} for it. Let us suppose that this conditional expectation is a linear function of z. Then, for example, referring to Ning et al. (2020), the estimator is given by solving We will write this asα CB and its limit as α † . Subsequently,θ CB is obtained by solving It can be seen that the expectation of the left-hand side of (32) becomes 0 p if E{(∂/∂θ) log f (y (h) | x (h) ; θ) | z} is a linear function of z. This means thatθ CB converges in probability to θ * . In this setting, let us derive a criterion like IPWIC in Section 3.3, for example, for the selection of x. It differs from the one in Section 3.3 in that it supposes that α † is not necessarily the true value. when the confounding variables are high dimensional is a timely topic in econometrics (e.g., Belloni et al. 2017, Chernozhukov et al. 2018, Athey et al. 2018, but even in this setting, it is obvious that, for an efficient estimation, we must select regression models for the explanatory variables and select confounding variables that directly affect the causal effects themselves when they are not so numerous. Extending the information criterion to deal with this problem is also a future challenge.