constraints.affine¶
Module: constraints.affine¶
Inheritance diagram for selectinf.constraints.affine:
This module contains the core code needed for post selection inference based on affine selection procedures as described in the papers Kac Rice, Spacings, covTest and post selection LASSO.
Classes¶
constraints¶
-
class
selectinf.constraints.affine.constraints(linear_part, offset, covariance=None, mean=None, rank=None)[source]¶ Bases:
objectThis class is the core object for affine selection procedures. It is meant to describe sets of the form \(C\) where
\[C = \left\{z: Az\leq b \right \}\]Its main purpose is to consider slices through \(C\) and the conditional distribution of a Gaussian \(N(\mu,\Sigma)\) restricted to such slices.
Notes
In this parameterization, the parameter self.mean corresponds to the reference measure that is being truncated. It is not the mean of the truncated Gaussian.
>>> positive = constraints(-np.identity(2), np.zeros(2)) >>> Y = np.array([3, 4.4]) >>> eta = np.array([1, 1], np.float) >>> list(positive.interval(eta, Y)) [4.62..., 10.17...] >>> positive.pivot(eta, Y) 5.187...-07 >>> list(positive.bounds(eta, Y)) [1.399..., 7.4..., inf, 1.414...] >>>
-
__init__(linear_part, offset, covariance=None, mean=None, rank=None)[source]¶ Create a new inequality.
- Parameters
linear_part : np.float((q,p))
The linear part, \(A\) of the affine constraint \(\{z:Az \leq b\}\).
offset: np.float(q)
The offset part, \(b\) of the affine constraint \(\{z:Az \leq b\}\).
covariance : np.float((p,p))
Covariance matrix of Gaussian distribution to be truncated. Defaults to np.identity(self.dim).
mean : np.float(p)
Mean vector of Gaussian distribution to be truncated. Defaults to np.zeros(self.dim).
rank : int
If not None, this should specify the rank of the covariance matrix. Defaults to self.dim.
-
conditional(linear_part, value, rank=None)[source]¶ Return an equivalent constraint after having conditioned on a linear equality.
Let the inequality constraints be specified by (A,b) and the equality constraints be specified by (C,d). We form equivalent inequality constraints by considering the residual
\[AY - E(AY|CY=d)\]- Parameters
linear_part : np.float((k,q))
Linear part of equality constraint, C above.
value : np.float(k)
Value of equality constraint, b above.
rank : int
If not None, this should specify the rank of linear_part. Defaults to min(k,q).
- Returns
conditional_con : constraints
Affine constraints having applied equality constraint.
-
bounds(direction_of_interest, Y)[source]¶ For a realization \(Y\) of the random variable \(N(\mu,\Sigma)\) truncated to \(C\) specified by self.constraints compute the slice of the inequality constraints in a given direction \(\eta\).
- Parameters
direction_of_interest: np.float
A direction \(\eta\) for which we may want to form selection intervals or a test.
Y : np.float
A realization of \(N(\mu,\Sigma)\) where \(\Sigma\) is self.covariance.
- Returns
L : np.float
Lower truncation bound.
Z : np.float
The observed \(\eta^TY\)
U : np.float
Upper truncation bound.
S : np.float
Standard deviation of \(\eta^TY\).
-
pivot(direction_of_interest, Y, null_value=None, alternative='greater')[source]¶ For a realization \(Y\) of the random variable \(N(\mu,\Sigma)\) truncated to \(C\) specified by self.constraints compute the slice of the inequality constraints in a given direction \(\eta\) and test whether \(\eta^T\mu\) is greater then 0, less than 0 or equal to 0.
- Parameters
direction_of_interest: np.float
A direction \(\eta\) for which we may want to form selection intervals or a test.
Y : np.float
A realization of \(N(0,\Sigma)\) where \(\Sigma\) is self.covariance.
alternative : [‘greater’, ‘less’, ‘twosided’]
What alternative to use.
- Returns
P : np.float
\(p\)-value of corresponding test.
Notes
All of the tests are based on the exact pivot \(F\) given by the truncated Gaussian distribution for the given direction \(\eta\). If the alternative is ‘greater’ then we return \(1-F\); if it is ‘less’ we return \(F\) and if it is ‘twosided’ we return \(2 \min(F,1-F)\).
-
interval(direction_of_interest, Y, alpha=0.05, UMAU=False)[source]¶ For a realization \(Y\) of the random variable \(N(\mu,\Sigma)\) truncated to \(C\) specified by self.constraints compute the slice of the inequality constraints in a given direction \(\eta\) and test whether \(\eta^T\mu\) is greater then 0, less than 0 or equal to 0.
- Parameters
direction_of_interest: np.float
A direction \(\eta\) for which we may want to form selection intervals or a test.
Y : np.float
A realization of \(N(0,\Sigma)\) where \(\Sigma\) is self.covariance.
alpha : float
What level of confidence?
UMAU : bool
Use the UMAU intervals?
- Returns
[U,L] : selection interval
-
covariance_factors(force=True)[source]¶ Factor self.covariance, finding a possibly non-square square-root.
- Parameters
force : bool
If True, force a recomputation of the covariance. If not, assumes that covariance has not changed.
-
estimate_mean(observed)[source]¶ Softmax estimator based on an observed data point.
Makes a whitened copy then returns softmax estimate.
TODO: what if self.mean != 0 before hand?
-
gaussian_hit_and_run¶
-
class
selectinf.constraints.affine.gaussian_hit_and_run(constraints, state, nstep=1)[source]¶ Bases:
selectinf.distributions.chain.reversible_markov_chain-
__init__(constraints, state, nstep=1)[source]¶ Initialize self. See help(type(self)) for accurate signature.
-
backward_step()¶ For a reversible chain, a backward_step is just a step.
-
forward_step()¶ For a reversible chain, a forward_step is just a step.
-
get_state()¶
-
next()¶
-
set_state(state)¶
-
property
state¶
-
Functions¶
-
selectinf.constraints.affine.gibbs_test(affine_con, Y, direction_of_interest, how_often=-1, ndraw=5000, burnin=2000, white=False, alternative='twosided', UMPU=True, sigma_known=False, alpha=0.05, pvalue=True, use_constraint_directions=False, use_random_directions=True, tilt=None, test_statistic=None, accept_reject_params=(100, 15, 2000), MLE_opts={'burnin': 1000, 'hessian_min': 1.0, 'how_often': 5, 'ndraw': 500, 'niter': 10, 'startMLE': None, 'step_size': 0.2, 'tol': 1e-06})[source]¶ A Monte Carlo significance test for a given function of con.mean.
- Parameters
affine_con : `selection.affine.constraints`_
Y : np.float
Point satisfying the constraint.
direction_of_interest: np.float
Which linear function of con.mean is of interest? (a.k.a. \(\eta\) in many of related papers)
how_often : int (optional)
How often should the sampler make a move along direction_of_interest? If negative, defaults to ndraw+burnin (so it will never be used).
ndraw : int (optional)
Defaults to 5000.
burnin : int (optional)
Defaults to 2000.
white : bool (optional)
Is con.covariance equal to identity?
alternative : str
One of [‘greater’, ‘less’, ‘twosided’]
UMPU : bool
Perform the UMPU test?
sigma_known : bool
Is \(\sigma\) assumed known?
alpha :
Level for UMPU test.
pvalue :
Return a pvalue or just the decision.
use_constraint_directions : bool (optional)
Use the directions formed by the constraints as in the Gibbs scheme?
use_random_directions : bool (optional)
Use additional random directions in the Gibbs scheme?
tilt : np.float (optional)
If not None, a direction to tilt. The likelihood ratio is effectively multiplied by \(y^TP\eta\) where \(\eta\) is tilt, and \(P\) is projection onto the rowspace of affine_con.
accept_reject_params : tuple
If not () should be a tuple (num_trial, min_accept, num_draw). In this case, we first try num_trial accept-reject samples, if at least min_accept of them succeed, we just draw num_draw accept_reject samples.
MLE_opts : {}
Arguments passed to one_parameter_MLE if tilt is not None.
- Returns
pvalue : float
P-value (using importance weights) for specified hypothesis test.
Z : np.float((ndraw, n))
Sample from the sphere intersect the constraints.
weights : np.float(ndraw)
Importance weights for the sample.
dfam : discrete_family
Discrete exponential family with above sufficient statistics Z and weights.
-
selectinf.constraints.affine.interval_constraints(support_directions, support_offsets, covariance, observed_data, direction_of_interest, tol=0.0001)[source]¶ Given an affine constraint \(\{z:Az \leq b \leq \}\) (elementwise) specified with \(A\) as support_directions and \(b\) as support_offset, a new direction of interest \(\eta\), and an observed_data is Gaussian vector \(Z \sim N(\mu,\Sigma)\) with covariance matrix \(\Sigma\), this function returns \(\eta^TZ\) as well as an interval bounding this value.
The interval constructed is such that the endpoints are independent of \(\eta^TZ\), hence the \(p\)-value of Kac Rice can be used to form an exact pivot.
- Parameters
support_directions : np.float
Matrix specifying constraint, \(A\).
support_offsets : np.float
Offset in constraint, \(b\).
covariance : np.float
Covariance matrix of observed_data.
observed_data : np.float
Observations.
direction_of_interest : np.float
Direction in which we’re interested for the contrast.
tol : float
Relative tolerance parameter for deciding sign of \(Az-b\).
- Returns
lower_bound : float
observed : float
upper_bound : float
sigma : float
-
selectinf.constraints.affine.sample_from_constraints(con, Y, direction_of_interest=None, how_often=-1, ndraw=1000, burnin=1000, white=False, use_constraint_directions=True, use_random_directions=True, accept_reject_params=())[source]¶ Use Gibbs sampler to simulate from con.
- Parameters
con : `selection.affine.constraints`_
Y : np.float
Point satisfying the constraint.
direction_of_interest : np.float (optional)
Which projection is of most interest?
how_often : int (optional)
How often should the sampler make a move along direction_of_interest? If negative, defaults to ndraw+burnin (so it will never be used).
ndraw : int (optional)
Defaults to 1000.
burnin : int (optional)
Defaults to 1000.
white : bool (optional)
Is con.covariance equal to identity?
use_constraint_directions : bool (optional)
Use the directions formed by the constraints as in the Gibbs scheme?
use_random_directions : bool (optional)
Use additional random directions in the Gibbs scheme?
accept_reject_params : tuple
If not () should be a tuple (num_trial, min_accept, num_draw). In this case, we first try num_trial accept-reject samples, if at least min_accept of them succeed, we just draw num_draw accept_reject samples.
- Returns
Z : np.float((ndraw, n))
Sample from the Gaussian distribution conditioned on the constraints.
-
selectinf.constraints.affine.sample_from_sphere(con, Y, direction_of_interest=None, how_often=-1, ndraw=1000, burnin=1000, use_constraint_directions=True, use_random_directions=True, white=False)[source]¶ Use Gibbs sampler to simulate from con intersected with (whitened) sphere of radius np.linalg.norm(Y). When con.covariance is not \(I\), it samples from the ellipse of constant Mahalanobis distance from con.mean.
- Parameters
con : `selection.affine.constraints`_
Y : np.float
Point satisfying the constraint.
direction_of_interest : np.float (optional)
Which projection is of most interest?
how_often : int (optional)
How often should the sampler make a move along direction_of_interest? If negative, defaults to ndraw+burnin (so it will never be used).
ndraw : int (optional)
Defaults to 1000.
burnin : int (optional)
Defaults to 1000.
white : bool (optional)
Is con.covariance equal to identity?
- Returns
Z : np.float((ndraw, n))
Sample from the sphere intersect the constraints.
weights : np.float(ndraw)
Importance weights for the sample.
-
selectinf.constraints.affine.selection_interval(support_directions, support_offsets, covariance, observed_data, direction_of_interest, tol=0.0001, alpha=0.05, UMAU=True)[source]¶ Given an affine in cone constraint \(\{z:Az+b \leq 0\}\) (elementwise) specified with \(A\) as support_directions and \(b\) as support_offset, a new direction of interest \(\eta\), and an observed_data is Gaussian vector \(Z \sim N(\mu,\Sigma)\) with covariance matrix \(\Sigma\), this function returns a confidence interval for \(\eta^T\mu\).
- Parameters
support_directions : np.float
Matrix specifying constraint, \(A\).
support_offset : np.float
Offset in constraint, \(b\).
covariance : np.float
Covariance matrix of observed_data.
observed_data : np.float
Observations.
direction_of_interest : np.float
Direction in which we’re interested for the contrast.
tol : float
Relative tolerance parameter for deciding sign of \(Az-b\).
UMAU : bool
Use the UMAU interval, or twosided pivot.
- Returns
selection_interval : (float, float)
-
selectinf.constraints.affine.stack(*cons)[source]¶ Combine constraints into a large constaint by intersection.
- Parameters
cons : [`selection.affine.constraints`_]
A sequence of constraints.
- Returns
intersection : `selection.affine.constraints`_
Notes
Resulting constraint will have mean 0 and covariance \(I\).