# A note on estimation in a simple probit model under dependency

We consider a probit model without covariates, but the latent Gaussian variables having compound symmetry covariance structure with a single parameter characterizing the common correlation. We study the parameter estimation problem under such one-parameter probit models. As a surprise, we demonstrate that the likelihood function does not yield consistent estimates for the correlation. We then formally prove the parameter's nonestimability by deriving a non-vanishing minimax lower bound. This counter-intuitive phenomenon provides an interesting insight that one bit information of the latent Gaussian variables is not sufficient to consistently recover their correlation. On the other hand, we further show that trinary data generated from the Gaussian variables can consistently estimate the correlation with parametric convergence rate. Hence we reveal a phase transition phenomenon regarding the discretization of latent Gaussian variables while preserving the estimability of the correlation.

## Authors

• 5 publications
• 57 publications
• ### Nonparametric Latent Tree Graphical Models: Inference, Estimation, and Structure Learning

Tree structured graphical models are powerful at expressing long range o...
01/16/2014 ∙ by Le Song, et al. ∙ 0

• ### Bayesian Sparse Covariance Structure Analysis for Correlated Count Data

In this paper, we propose a Bayesian Graphical LASSO for correlated coun...
06/05/2020 ∙ by Sho Ichigozaki, et al. ∙ 0

• ### Estimation in Tensor Ising Models

The p-tensor Ising model is a one-parameter discrete exponential family ...
08/29/2020 ∙ by Somabha Mukherjee, et al. ∙ 0

• ### Solving Non-identifiable Latent Feature Models

Latent feature models (LFM)s are widely employed for extracting latent s...
09/11/2018 ∙ by Ryota Suzuki, et al. ∙ 0

• ### Variational inference for sparse network reconstruction from count data

In multivariate statistics, the question of finding direct interactions ...
06/08/2018 ∙ by Julien Chiquet, et al. ∙ 0

• ### Synthesis of Gaussian Trees with Correlation Sign Ambiguity: An Information Theoretic Approach

In latent Gaussian trees the pairwise correlation signs between the vari...
01/24/2016 ∙ by Ali Moharrer, et al. ∙ 0

• ### A Note on Implementing a Special Case of the LEAR Covariance Model in Standard Software

Repeated measures analyses require proper choice of the correlation mode...
07/26/2017 ∙ by Sean L. Simpson, et al. ∙ 0

##### This week in AI

Get the week's most popular data science and artificial intelligence research sent straight to your inbox every Saturday.

## 1 Problem statement

We consider a simple one-parameter probit model under dependency. The observed binary data

is generated by thresholding a multivariate Gaussian vector with compound symmetry covariance structure:

 Ai=1(Xi>0),  i=1,…,n,(X1,…,Xn)∼N(0,Σ), (1) Σii=1,Σij=a∗,  i=1,…,n,j=1,…,n,i≠j.

The single parameter characterizes the dependency in the exchangeable sequence . The covariance matrix can be considered as a special case of the spiked covariance form proposed in Johnstone & Lu (2004, 2009)

for studying high-dimensional principal component analysis. The central question we focus on in this paper is

As , can be consistently estimated from the sequence ?

In contrast to usual probit models where and the mean of is parameterized by a function of covariates, model (1) puts emphasis on the dependency structure of the observations, and is possibly one of the simplest models along this line.

The model can be well motivated from network data analysis. For notational convenience, we rewrite the observations by a matrix and the latent Gaussian variables by . The Bernoulli entry represents whether there exists an edge between nodes and . We consider undirected networks where for simplicity. The covariance matrix belongs to . Due to rich structures exhibited in different types of networked systems (Newman, 2003), there has been an extensive literature on network modeling including Erdös–Rényi random graph model (Erdös & Rényi, 1959), exponential random graph model (Robins et al., 2007), stochastic blockmodel (Holland et al., 1983), and latent space model (Hoff et al., 2002), among others. As an alternative, model (1) assumes the edges between nodes are generated from underlying Gaussian variables, and the covariance matrix captures the possible dependency among different edges. When , (1) is an example of Erdös–Rényi random graph. We should emphasize that model (1) is not sophisticated enough for fitting most real network data. One possible generalization is

 Aij=1(Xij>0),i=1,…,m,j=1,…,m, X=(Xij)m×m∼N(Θ,Σ).

The mean parameter is introduced to incorporate heterogeneity across different edges. It can be assumed of low-rank based on the hypothesis that the generation mechanism of edges is driven by a few node-specific factors. This is in the same spirit of both stochastic block model and latent space model. Furthermore, more general structure can be imposed for . For instance, letting be the covariance between and , then one possible generalization is

 Σ(i,j),(k,l)=⎧⎨⎩1{i,j}∩{k,l}={i,j}∪{k,l},b∗{i,j}∩{k,l}=∅,a∗otherwise, (2)

where is the empty set. The formulation (2) automatically enforces the symmetry constraint which is required for modeling undirected networks. In this case, the dependencies among edges that share common nodes and those which do not are characterized by two different parameters. To fix idea, we will focus on the simplified model (1). Nevertheless, our analysis will shed lights on general probit modeling of networks. We start by investigating the likelihood approach for estimating and indicating its infeasibility. We then formally prove the nonestimability of under model (1) and provide a solution to make the estimation possible. Finally we discuss some implications of our results regarding binary data modeling.

## 2 Likelihood methods

Since model (1) takes a simple parametric form, our first attempt is to check if the maximum likelihood estimate is consistent. Before that we give an alternative formulation of the model that will be useful in later discussions.

###### Lemma 1.

Model (1) can be reformulated as

 Ai=1((1−a∗)1/2Yi+(a∗)1/2Y>0),  i=1,…,n, (3)

where are independently and identically distributed as .

Let

be the cumulative distribution function of

, , and

 Fn(z)=¯AlogΦ((a)1/2(1−a)−1/2z)+(1−¯A)logΦ(−(a)1/2(1−a)−1/2z).

According to Lemma 1, we can write down the likelihood function,

 Ln(a;{Ai}) = ∫∞−∞[Φ((a)1/2(1−a)−1/2z)]∑ni=1Ai[Φ(−(a)1/2(1−a)−1/2z)]n−∑ni=1Aiexp(−z2/2)dz = ∫∞−∞exp(nFn(z)−z2/2)dz.

Ideally, we hope a properly normalized converges to a deterministic function whose maximizer is . But the following result says this is not the case.

###### Proposition 1.

Under the model formulation (3), for any

, with probability 1,

 n−1logLn(a;{Ai}) → Φ((a∗)1/2(1−a∗)−1/2Y)logΦ((a∗)1/2(1−a∗)−1/2Y)+ Φ(−(a∗)1/2(1−a∗)−1/2Y)logΦ(−(a∗)1/2(1−a∗)−1/2Y),

as .

Proposition 1 is essentially a first order Laplacian approximation. It shows that the normalized log-likelihood function in (1) does not converge to a deterministic function, and more interestingly the limiting function is invariant of for . Figure 1 illustrates the normalized log-likelihood function curves when and the sample size . It is clear that the functions are pretty flat in most of . Interestingly, the functions have a sharp transition at around . This can be justified by a straightforward calculation showing . Hence the limit at is deterministic and does not even depend on .

While the normalized likelihood function converges (in first order) to an unfavorable object, it does not necessarily mean that likelihood based estimates are not consistent, because the higher order terms in the limiting function may contain useful information for the parameter of interest. Unfortunately, this is not true either, as demonstrated below by a second order analysis of the likelihood function.

###### Proposition 2.

Under the model formulation (3), for any , as ,

 Ln(a;{Ai})exp(−n(¯Alog¯A+(1−¯A)log(1−¯A))) = [2π(1−a)a−1Φ((a∗1−a∗)1/2Y)Φ(−(a∗1−a∗)1/2Y)]12ψ((a∗1−a∗)1/2Y)exp(a∗(a−1)Y22(1−a∗)a)n−1/2+Op(n−1),

where

is the probability density function of

.

Proposition 2 can be considered as a second order Laplacian approximation result. Proposition 1 shows that the first order term of is exponentially small. According to Proposition 2, after is scaled by a fully data-dependent and exponentially small term , the dominating term is of order . Moreover, the number comes into play in the precise constant of the second order term of , and the constant attains maximum at . If we replace with its expectation, then would equal . However, since the second order term is a random function of , we may not expect maximum likelihood estimate to be consistent for .

## 3 Nonestimability of a∗

The preceding discussions in Section 2 indicate that the likelihood approach may not yield consistent estimates for . Since model (1

) is a simple parametric model, it may further imply no consistent estimates exist. Indeed, this is formally proved in the theorem below.

###### Theorem 1.

Under model (1), given any strictly increasing function with and any ,

 liminfn→∞infTmaxa∗∈{a1,a2}Ea∗[R(|T({Ai}ni=1)−a∗|)]>0,

where is any measurable function and the expectation is taken over .

###### Remark 1.

Theorem 1 reveals that no consistent estimates of exist even when the parameter space only consists of two elements. At first glance, this result seems counter-intuitive. The model (1) has a simple structure with a single parameter, but there is no way to reliably estimate the parameter from infinite components of the data sequence. On the other hand, if the observation is the hidden sequence instead of , -consistent estimates can be easily constructed. For instance, denote . Lemma 1 implies that are independently and identically distributed as . Hence is -consistent for .

###### Remark 2.

One might argue that it is the strong dependency in causing the issue. If the dependency can be somehow weakened, will become estimable? For example, we can consider with a constant and , and ask if there is any estimate such that . With a similar proof as the one for Theorem 1, it is possible to show

 liminfn→∞infTmaxa∗0∈{a1,a2}Ea∗(ρ−1n|T({Ai}ni=1)−a∗|)>0.

Hence remains nonestimable.

## 4 A solution

The arguments in the last two paragraphs suggest that the nonestimability of in model (1) is due to the thresholding operation in the binary data generating process. Specifically, is estimable from , but the thresholding on to produce results in loss of too much information, to be able to consistently estimate . In a nutshell, one bit information of each is not sufficient to recover the common correlation among them.

Interestingly, we demonstrate below that in fact a little more than one bit of information is enough to obtain -consistent estimates. Given two constants , denote . Suppose

 Ai=(A(1)i,A(2)i,A(3)i), A(j)i=1(Xi∈Ij), j=1,2,3, i=1,…,n.

That is, for every , instead of observing the sign of , we know which one of the three intervals belongs to. When , it reduces to the model (1). The theorem below gives one -consistent estimate.

###### Theorem 2.

For all , consider the estimate

 ^an=1−[t1−t2Φ−1(¯A(1))+Φ−1(¯A(3))]21(0<¯A(1)<1,0<¯A(3)<1), (5)

where , and is the inverse function of . Then as , .

Figure 2 empirically verifies the consistency result of in Theorem 2.

###### Remark 3.

Following the preceding discussion, we can further consider a general setting where there exist () consecutive intervals and represents which interval falls into. It is not hard to show that -consistent estimate for can be constructed in a similar way as in (5). At a high level, we may consider the observed sequence as a discretized version of . Theorems 1 and 2 together characterize a phase transition regarding the estimability of . That is, is estimable if and only if .

## 5 Discussion

We have revealed an interesting phenomenon regarding the estimability of a single parameter in a simple probit model. Several important directions are left open. For instance, suppose the latent vector is modeled by the Gaussian copula family (Klaassen et al., 1997; Tsukahara, 2005) with the parameter having compound symmetry covariance structure as in (1). What can we say about the estimability of ? As another example, consider the covariance matrix in (1) is replaced by the more general one in (2) in network modeling. It is clear that both and are nonestimable from the one bit information of . The question is how much more information is needed to consistently estimate them. Would the phase transition phenomenon we discussed in Remark 3 continue to hold?

Our results also have a few insightful implications. For example, modeling a dependent exchangeable binary sequence is subtle. We have shown that it is even impossible to do inference under a simple one-parameter model. Regarding binary network modeling, it might not be desirable to assume dependency among all the edges. Furthermore, converting a weighted network to a binary one may lose substantial information for the sake of parameter estimation.

## Acknowledgement

The authors are grateful to Professor Anthony C. Davison and Professor Zhiliang Ying for their insightful comments which greatly improved the scope and presentation of this paper.

## Appendix A Notations

Recall that are the cumulative distribution function and probability density function of a standard normal respectively, and is the inverse function of . We will use the following notations extensively,

 ¯A = n−1n∑i=1Ai,  z∗=(1−a)1/2a−1/2Φ−1(¯A), Fn(z) = ¯AlogΦ(a1/2(1−a)−1/2z)+(1−¯A)logΦ(−a1/2(1−a)−1/2z).

We also use to mean and are of the same order. The next section collects a few useful lemmas that will be applied in the later proofs.

## Appendix B Useful lemmas

###### Lemma 2.

The following holds

(i) is strictly log-concave
(ii)
(iii)
(iv)
(v) , as .

###### Proof.

Part (i) can be easily verified. Parts (ii) and (iii) can be shown by using the tail probability , as . Part (iv) is due to the fact that the second derivative of is negative for any and goes to only when . Part (v) is taken from (2.6) in Wong (2001) with minor change. ∎

###### Lemma 3.

Consider model (1), as ,

 pr(¯A=0)=pr(¯A=1)∝na∗−1a∗(logn)1−2a∗2a∗.
###### Proof.

The proof is a direct adaption from the calculations of the integral (2.1) in Wong (2001). We hence do not repeat the arguments. ∎

###### Lemma 4.

Consider model (1), as ,

 E[(Φ−1(¯A))21(0<¯A<1)]=O(1),E[(−log[(1−¯A)¯A])1(0<¯A<1)]=O(1).
###### Proof.

From Lemma 2 part (v), it is straightforward to confirm that

 limx→0(Φ−1(x))2−log(x(1−x))=limx→1(Φ−1(x))2−log(x(1−x))=2.

So there exists such that for . Hence it suffices to show is bounded. From the model formulation (3), it is easy to see that

 E[−log¯A1(0<¯A<1)]=E[−log(1−¯A)1(0<¯A<1)].

We thus focus on bounding . For notational simplicity, we denote

 p(Y)=E(Ai∣Y)=Φ((a∗)1/2(1−a∗)−1/2Y).

First,

 E[−log¯A1(0<¯A<1)] ≤ E[−log¯A1(1/n≤¯A

Also, for large ,

 pr(¯A

where we have used Hoeffding’s inequality in , and the inequality for large in . Hence we obtain . The proof will be completed if we can show . This is true because when is large enough. ∎

## Appendix C Proof of Proposition 1

###### Proof.

It is straightforward to verify that when . Hence we can readily have the upper bound, if

 logLn(a;{Ai})n=1nlog∫∞−∞exp(nFn(z)−z2/2)dz≤Fn(z∗)+12nlog(2π). (6)

Regarding the lower bound, if there exists an absolute constant such that

 1nlog∫∞−∞exp(nFn(z)−z2/2)dz (a)= 1nlog∫∞−∞exp(nFn(z∗)+nF′′n(~z)(z−z∗)2/2−z2/2)dz (b)≥ Fn(z∗)+1nlog∫∞−∞exp(−C1na2(1−a)(z−z∗)2−z22)dz ≥ Fn(z∗)+log(2π)2n−(z∗)22n−log(C1na1−a+1)2n,

where is by a Taylor expansion and is due to Lemma 2 part (ii). Lemma 1 from the main text implies that , almost surely. Hence, almost surely as ,

 z∗→(a∗(1−a))1/2(a(1−a∗))−1/2Y,1(0<¯A<1)→1, Fn(z∗)→Φ((a∗)1/2Y(1−a∗)1/2)logΦ((a∗)1/2Y(1−a∗)1/2)+Φ(−(a∗)1/2Y(1−a∗)1/2)logΦ(−(a∗)1/2Y(1−a∗)1/2).

These results combined with the upper and lower bounds in (6) and (C) completes the proof. ∎

## Appendix D Proof of Proposition 2

###### Proof.

We first restrict our analysis to the case . We have

 Ln(a;{Ai})exp(−n(¯Alog¯A+(1−¯A)log(1−¯A))) = ∫z∗−∞exp(n(Fn(z)−Fn(z∗))−z2/2)dz+∫∞z∗exp(n(Fn(z)−Fn(z∗))−z2/2)dz.

We focus on the first integral above. Lemma 2 part (i) says that is strictly concave. Hence is increasing in . By a change of variable we obtain

 ∫z∗−∞exp(n(Fn(z)−Fn(z∗))−z2/2)dz=∫∞0exp(−n~z)exp(−z2/2)(F′n(z))−1g(~z)d~z.

Denote . Then yielding

 ∫z∗−∞exp(n(Fn(z)−Fn(z∗))−z2/2)dz (8) = h(0)(πn)1/2+∫10∫∞0exp(−n~z)h′(t~z)(~z)1/2d~zdt.

Based on the Taylor expansions

 ~z=Fn(z∗)−Fn(z)=−12F′′n(z1)(z−z∗)2,F′n(z)=F′′n(z2)(z−z∗) (9)

with , the following holds

 h(0) = lim~z→0exp(−z2/2)(F′n(z))−1(~z)1/2 = limz↗z∗exp(−z2/2)[−12F′′n(z1)(z−z∗)2]1/2F′′n(z2)(z−z∗)=exp(−(z∗)2/2)(−2F′′n(z∗))1/2.

Next we bound the second term on the right-hand side of (8). With a few calculations we get

 h′(~z)(~z)1/2 = exp(−z2/2)z(Fn(z∗)−Fn(z))(F′n(z))2+ (10) exp(−z2/2)F′′n(z)(Fn(z∗)−Fn(z))+12(F′n(z))2(F′n(z))3.

According to Lemma 2 part (ii) and (iv), there exist absolute constants such that for any

 −a1−aC1≤F′′n(z) (11) = a1−a{¯A[logΦ((a)1/2z/(1−a)1/2)]′′+(1−¯A)[logΦ(−(a)1/2z/(1−a)1/2)]′′} ≤ −a1−aC2min(¯A,1−¯A).

Then using the expansions we had about and in (9), we can obtain

 supzexp(−z2/2)|z(Fn(z∗)−Fn(z))|(F′n(z))2=supzexp(−z2/2)|zF′′n(z1)|2(F′′n(z2))2≤C3min(¯A,1−¯A)2, (12)

where is a positive constant only depending on . To bound the other term in (10), we need take the Taylor expansions to a higher order,

 F′n(z)=F′′n(z∗)(z−z∗)+12F′′′n(z4)(z−z∗)2,  z4∈[z,z∗], (13) Fn(z∗)−Fn(z)=−12F′′n(z∗)(z−z∗)2−16F′′′n(z3)(z−z∗)3,  z3∈[z,z∗]. (14)

Plugging the above expansions together with in the numerator of the second term in (10) and use in the denominator, it is not hard to obtain

 supz∣∣exp(−z2/2)F′′n(z)(Fn(z∗)−Fn(z))+12(F′n(z))2(F′n(z))3∣∣ (a)≤ supzexp(−z2/2)(C4+C5(|z|+|z∗|))|F′′n(z2)|3(b)≤C6+C7|z∗|min(¯A,1−¯A)3,

where the positive constants only depend on . To derive we have used Lemma 2 part (ii) (iii), and is due to (11). Combining the above upper bound with (