​
In[]:=
$Version
Out[]=
12.1.0 for Mac OS X x86 (64-bit) (March 14, 2020)
In[]:=
$OperatingSystem
Out[]=
MacOSX
The idea of this code is to find a correlation matrix such that the partial correlation between a (2nd feature) and y (what you are trying to predict) has the potential to be 0 even if the zero-order correlation is non-zero. The θ parameter controls whether the partial correlation is negative (θ<0) positive (θ>0) or zero (θ=0).
In[]:=
generateCorrelationMatrix[rxa_,rxy_,θ_/;-1≤θ≤1]:=Module[{product=rxa*rxy,ray,matrix},ray=Piecewise[{{Rescale[θ,{0,1},{product,1}],0≤θ≤1},{Rescale[θ,{-1,0},{0,product}],-1≤θ<0}}];​​matrix={{1,rxa,rxy},{rxa,1,ray},{rxy,ray,1}};​​If[Not[PositiveDefiniteMatrixQ[matrix]],$Failed,matrix]​​]​​
Let' s generate a copula with θ = 0.
In[]:=
copula=CopulaDistribution[{"Multinormal",generateCorrelationMatrix[0.3,0.7,0.]},{NormalDistribution[],BernoulliDistribution[0.75],BernoulliDistribution[0.8]}]
Out[]=
CopulaDistribution[{Multinormal,{{1,0.3,0.7},{0.3,1,0.21},{0.7,0.21,1}}},{NormalDistribution[0,1],BernoulliDistribution[0.75],BernoulliDistribution[0.8]}]
Seed our random number generator:
In[]:=
SeedRandom[2]
generate some data
In[]:=
data=RandomVariate[copula,20]
Out[]=
{{0.623366,1,1},{-0.272176,1,1},{-1.52651,0,0},{-0.744493,1,1},{0.400523,1,1},{0.404079,1,1},{-0.724666,1,1},{-0.750685,0,1},{-0.782566,0,0},{-0.687976,1,0},{0.432732,1,1},{-0.159116,1,1},{-0.00144431,1,0},{-1.13231,0,1},{-1.34211,1,0},{1.274,0,1},{0.368973,0,1},{-0.86363,0,1},{1.58101,1,1},{0.146163,1,1}}
run logistic regression
In[]:=
fail=LogitModelFit[data,{1,x,a},{x,a},NominalVariables{a}]
General
:Exp[-979888.] is too small to represent as a normalized machine number; precision may be lost.
General
:Exp[-979856.] is too small to represent as a normalized machine number; precision may be lost.
General
:Exp[-979856.] is too small to represent as a normalized machine number; precision may be lost.
General
:Further output of General::munfl will be suppressed during this calculation.
LogitModelFit
:Non-real values were encountered during the IRLS process. The specified model may not be appropriate for the data.
Out[]=
LogitModelFit[{{0.623366,1,1},{-0.272176,1,1},{-1.52651,0,0},{-0.744493,1,1},{0.400523,1,1},{0.404079,1,1},{-0.724666,1,1},{-0.750685,0,1},{-0.782566,0,0},{-0.687976,1,0},{0.432732,1,1},{-0.159116,1,1},{-0.00144431,1,0},{-1.13231,0,1},{-1.34211,1,0},{1.274,0,1},{0.368973,0,1},{-0.86363,0,1},{1.58101,1,1},{0.146163,1,1}},{1,x,a},{x,a},NominalVariables{a}]
reseed and rerun
In[]:=
RandomSeed[9];​​jiggledData=RandomSample[data]
Out[]=
{{0.368973,0,1},{-0.782566,0,0},{1.274,0,1},{1.58101,1,1},{0.400523,1,1},{0.623366,1,1},{-0.687976,1,0},{-0.86363,0,1},{-1.52651,0,0},{-0.159116,1,1},{-0.750685,0,1},{0.146163,1,1},{-1.13231,0,1},{0.432732,1,1},{-0.00144431,1,0},{0.404079,1,1},{-0.272176,1,1},{-0.724666,1,1},{-1.34211,1,0},{-0.744493,1,1}}
In[]:=
success=LogitModelFit[jiggledData,{1,x,a},{x,a},NominalVariables{a}]
Out[]=
FittedModel
1
1+
-1.95656-1-1


In[]:=
success["ParameterTable"]
Out[]=
Estimate
Standard Error
z-Statistic
P-Value
1
1.95656
1.02401
1.91068
0.0560452
x
2.54915
1.39173
1.83164
0.0670044
a[0]
1.07365
1.45929
0.735733
0.461893