Skip to content

Commit

Permalink
styling
Browse files Browse the repository at this point in the history
  • Loading branch information
@fabian-sp
fabian-sp committed Nov 29, 2024
1 parent 7915562 commit aacbf98
Showing 1 changed file with 107 additions and 47 deletions.
154 changes: 107 additions & 47 deletions src/gglasso/problem.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -13,7 +13,6 @@

from .helper.model_selection import grid_search, single_grid_search, K_single_grid, ebic


assert_tol = 1e-5

class glasso_problem:
Expand Down Expand Up @@ -102,9 +101,14 @@ def __init__(self, S, N, reg = "GGL", reg_params = None, latent = False, G = Non


# create an instance of GGLassoEstimator (before scaling S!)
self.solution = GGLassoEstimator(S = self.S.copy(), N = self.N, p = self.p, K = self.K,\
multiple = self.multiple, latent = self.latent, conforming = self.conforming)

self.solution = GGLassoEstimator(S = self.S.copy(),
N = self.N,
p = self.p,
K = self.K,
multiple = self.multiple,
latent = self.latent,
conforming = self.conforming
)

# scale S by diagonal
if self.do_scaling:
Expand Down Expand Up @@ -141,8 +145,7 @@ def _derive_problem_formulation(self):
self.multiple = False

SINGLE_FLOAT_INT_TYPES = (int,float,np.int16,np.int32,np.int64,np.float16,np.float32,np.float64)



if type(self.S) == np.ndarray:

if len(self.S.shape) == 3:
Expand Down Expand Up @@ -425,29 +428,56 @@ def solve(self, Omega_0 = None, solver_params = dict(), tol = 1e-8, rtol = 1e-7,

if not self.multiple:
if self.latent:
sol, info = ADMM_SGL(S = self.S, lambda1 = self.reg_params['lambda1'], Omega_0 = self.Omega_0, \
tol = self.tol , rtol = self.rtol, latent = self.latent, mu1 = self.reg_params['mu1'],
lambda1_mask = self.reg_params.get('lambda1_mask'), **self.solver_params)
sol, info = ADMM_SGL(S=self.S,
lambda1=self.reg_params['lambda1'],
Omega_0=self.Omega_0,
tol=self.tol,
rtol=self.rtol,
latent=self.latent,
mu1=self.reg_params['mu1'],
lambda1_mask=self.reg_params.get('lambda1_mask'),
**self.solver_params
)

else:
sol = block_SGL(S = self.S, lambda1 = self.reg_params['lambda1'], Omega_0 = self.Omega_0, \
tol = self.tol, rtol = self.tol,
lambda1_mask = self.reg_params.get('lambda1_mask'), **self.solver_params)
sol = block_SGL(S=self.S,
lambda1=self.reg_params['lambda1'],
Omega_0=self.Omega_0,
tol=self.tol,
rtol=self.tol,
lambda1_mask=self.reg_params.get('lambda1_mask'),
**self.solver_params
)
info = {}


elif self.conforming:
sol, info = ADMM_MGL(S = self.S, lambda1 = self.reg_params['lambda1'], lambda2 = self.reg_params['lambda2'], reg = self.reg,\
Omega_0 = self.Omega_0, latent = self.latent, mu1 = self.reg_params['mu1'],\
tol = self.tol, rtol = self.rtol, **self.solver_params)


sol, info = ADMM_MGL(S=self.S,
lambda1=self.reg_params['lambda1'],
lambda2=self.reg_params['lambda2'],
reg=self.reg,
Omega_0=self.Omega_0,
latent=self.latent,
mu1=self.reg_params['mu1'],
tol=self.tol,
rtol=self.rtol,
**self.solver_params
)

else:
sol, info = ext_ADMM_MGL(S = self.S, lambda1 = self.reg_params['lambda1'], lambda2 = self.reg_params['lambda2'], reg = self.reg,\
Omega_0 = self.Omega_0, G = self.G, tol = self.tol, rtol = self.rtol,\
latent = self.latent, mu1 = self.reg_params['mu1'], **self.solver_params)


sol, info = ext_ADMM_MGL(S=self.S,
lambda1=self.reg_params['lambda1'],
lambda2=self.reg_params['lambda2'],
reg=self.reg,
Omega_0=self.Omega_0,
G=self.G,
tol=self.tol,
rtol=self.rtol,
latent=self.latent,
mu1=self.reg_params['mu1'],
**self.solver_params
)

# rescale
if self.do_scaling:
sol['Theta'] = self._rescale_to_covariances(sol['Theta'], self._scale)
Expand All @@ -457,9 +487,9 @@ def solve(self, Omega_0 = None, solver_params = dict(), tol = 1e-8, rtol = 1e-7,

# set the computed solution
if self.latent:
self.solution._set_solution(Theta = sol['Theta'], L = sol['L'])
self.solution._set_solution(Theta=sol['Theta'], L=sol['L'])
else:
self.solution._set_solution(Theta = sol['Theta'])
self.solution._set_solution(Theta=sol['Theta'])

self.solver_info = info.copy()
return
Expand All @@ -472,13 +502,13 @@ def solve(self, Omega_0 = None, solver_params = dict(), tol = 1e-8, rtol = 1e-7,
def _default_modelselect_params(self):

params = dict()
params['lambda1_range'] = np.logspace(0,-3,10)
params['lambda1_range'] = np.logspace(0, -3, 10)

if self.multiple:
params['lambda2_range'] = np.logspace(-1,-4,5)
params['lambda2_range'] = np.logspace(-1, -4, 5)

if self.latent:
params['mu1_range'] = np.logspace(2,-1,10)
params['mu1_range'] = np.logspace(2, -1, 10)
else:
params['mu1_range'] = None

Expand Down Expand Up @@ -510,7 +540,7 @@ def set_modelselect_params(self, modelselect_params = None):
warnings.warn("No grid for model selection is specified and thus default (or previous) values are used. A grid can be specified with the argument modelselect_params.")

else:
assert type(modelselect_params) == dict
assert isinstance(modelselect_params, dict)

# update with given input
# update with empty dict does not change the dictionary
Expand Down Expand Up @@ -573,10 +603,19 @@ def model_selection(self, modelselect_params = None, method = 'eBIC', gamma = 0.
# SINGLE GL --> GRID SEARCH lambda1/mu
###############################
if not self.multiple:
sol, self._all_theta, self._all_lowrank, stats = single_grid_search(S = self.S, lambda_range = self.modelselect_params['lambda1_range'], N = self.N, \
method = method, gamma = gamma, latent = self.latent, mu_range = self.modelselect_params['mu1_range'],
use_block = True, store_all = store_all, tol = tol, rtol = rtol,
lambda1_mask = self.modelselect_params['lambda1_mask'])
sol, self._all_theta, self._all_lowrank, stats = single_grid_search(S=self.S,
lambda_range=self.modelselect_params['lambda1_range'],
N=self.N,
method=method,
gamma=gamma,
latent=self.latent,
mu_range=self.modelselect_params['mu1_range'],
use_block=True,
store_all=store_all,
tol=tol,
rtol=rtol,
lambda1_mask=self.modelselect_params['lambda1_mask']
)

# update the regularization parameters to the best grid point
self.set_reg_params(stats['BEST'])
Expand All @@ -592,9 +631,18 @@ def model_selection(self, modelselect_params = None, method = 'eBIC', gamma = 0.
# LATENT VARIABLES --> FIRST STAGE lambda1/mu1 for each instance
###############################
if self.latent:
self._est_uniform, self._est_indv, stage1_statistics = K_single_grid(S = self.S, lambda_range = self.modelselect_params['lambda1_range'], N = self.N, method = method,\
gamma = gamma, latent = self.latent, mu_range = self.modelselect_params['mu1_range'],
use_block = True, store_all = store_all, tol = tol, rtol= rtol)
self._est_uniform, self._est_indv, stage1_statistics = K_single_grid(S=self.S,
lambda_range=self.modelselect_params['lambda1_range'],
N=self.N,
method=method,
gamma=gamma,
latent=self.latent,
mu_range=self.modelselect_params['mu1_range'],
use_block=True,
store_all=store_all,
tol=tol,
rtol=rtol
)

ix_mu = stage1_statistics['ix_mu']

Expand All @@ -608,10 +656,24 @@ def model_selection(self, modelselect_params = None, method = 'eBIC', gamma = 0.
# SECOND STAGE --> GRID SEARCH lambda1/lambda2
###############################

stats, best_ix, sol = grid_search(solver, S = self.S, N = self.N, p = self.p, reg = self.reg, l1 = self.modelselect_params['lambda1_range'], \
l2 = self.modelselect_params['lambda2_range'], w2 = None, method= method, gamma = gamma, \
G = self.G, latent = self.latent, mu_range = self.modelselect_params['mu1_range'], ix_mu = ix_mu, \
tol = tol, rtol = rtol, verbose = False)
stats, best_ix, sol = grid_search(solver,
S=self.S,
N=self.N,
p=self.p,
reg=self.reg,
l1=self.modelselect_params['lambda1_range'],
l2=self.modelselect_params['lambda2_range'],
w2=None,
method=method,
gamma=gamma,
G=self.G,
latent=self.latent,
mu_range=self.modelselect_params['mu1_range'],
ix_mu=ix_mu,
tol=tol,
rtol=rtol,
verbose=False
)

# update the lambda1/lambda2 regularization parameters to the best grid point
self.set_reg_params(stats['BEST'])
Expand All @@ -633,14 +695,12 @@ def model_selection(self, modelselect_params = None, method = 'eBIC', gamma = 0.

# set the computed solution
if self.latent:
self.solution._set_solution(Theta = sol['Theta'], L = sol['L'])
self.solution._set_solution(Theta=sol['Theta'], L=sol['L'])
else:
self.solution._set_solution(Theta = sol['Theta'])


self.solution._set_solution(Theta=sol['Theta'])

self.modelselect_stats = stats.copy()


return

#%%
Expand Down Expand Up @@ -711,18 +771,18 @@ def calc_ebic(self, gamma = 0.5):
"""
calculates the eBIC for a given value of :math:`\\gamma`. Note that this can differ from eBIC values in model selection because of the scaling.
"""
self.ebic_ = ebic(self.sample_covariance_, self.precision_, self.n_samples, gamma = gamma)
self.ebic_ = ebic(self.sample_covariance_, self.precision_, self.n_samples, gamma=gamma)
return self.ebic_

def calc_adjacency(self, t = 1e-8):

if self.conforming:
self.adjacency_ = adjacency_matrix(S = self.precision_, t = t)
self.adjacency_ = adjacency_matrix(S=self.precision_, t=t)

else:
self.adjacency_ = dict()
for k in range(self.K):
self.adjacency_[k] = adjacency_matrix(S = self.precision_[k], t = t)
self.adjacency_[k] = adjacency_matrix(S=self.precision_[k], t=t)

return

Expand Down

0 comments on commit aacbf98

Please sign in to comment.