Skip to content

Commit

Permalink
[pre-commit.ci] auto fixes from pre-commit.com hooks
Browse files Browse the repository at this point in the history
for more information, see https://pre-commit.ci
  • Loading branch information
pre-commit-ci[bot] committed Apr 8, 2024
1 parent 08f4583 commit e89265f
Show file tree
Hide file tree
Showing 35 changed files with 101 additions and 68 deletions.
93 changes: 47 additions & 46 deletions development/hessian_approximation/compute_hessian.py
@@ -1,6 +1,7 @@
"""This module provides the analytical solution for computing the hessian matrix of our
loglikelihood function
"""

import numpy as np
from scipy.stats import norm

Expand Down Expand Up @@ -29,9 +30,9 @@ def compute_hessian(x0, X1, X0, Z1, Z0, Y1, Y0):

# aux_params
nu1 = (Y1 - np.dot(beta1, X1.T)) / sd1
lambda1 = (np.dot(gamma, Z1.T) - rho1v * nu1) / (np.sqrt(1 - rho1v ** 2))
lambda1 = (np.dot(gamma, Z1.T) - rho1v * nu1) / (np.sqrt(1 - rho1v**2))
nu0 = (Y0 - np.dot(beta0, X0.T)) / sd0
lambda0 = (np.dot(gamma, Z0.T) - rho0v * nu0) / (np.sqrt(1 - rho0v ** 2))
lambda0 = (np.dot(gamma, Z0.T) - rho0v * nu0) / (np.sqrt(1 - rho0v**2))

eta1 = (
-lambda1 * norm.pdf(lambda1) * norm.cdf(lambda1) - norm.pdf(lambda1) ** 2
Expand Down Expand Up @@ -113,21 +114,21 @@ def calc_hess_beta1(
"""

# define some auxiliary variables
rho_aux1 = lambda1 * rho1v / (1 - rho1v ** 2) - nu1 / (1 - rho1v ** 2) ** 0.5
rho_aux2 = rho1v ** 2 / ((1 - rho1v ** 2) ** (3 / 2)) + 1 / (1 - rho1v ** 2) ** 0.5
sd_aux1 = rho1v ** 2 / (1 - rho1v ** 2)
sd_aux2 = rho1v / np.sqrt(1 - rho1v ** 2)
rho_aux1 = lambda1 * rho1v / (1 - rho1v**2) - nu1 / (1 - rho1v**2) ** 0.5
rho_aux2 = rho1v**2 / ((1 - rho1v**2) ** (3 / 2)) + 1 / (1 - rho1v**2) ** 0.5
sd_aux1 = rho1v**2 / (1 - rho1v**2)
sd_aux2 = rho1v / np.sqrt(1 - rho1v**2)
# derivation wrt beta1
der_b1_beta1 = -(
X1X1 * (rho1v ** 2 / (1 - rho1v ** 2)) * 1 / sd1 ** 2 - X1.T @ X1 / sd1 ** 2
X1X1 * (rho1v**2 / (1 - rho1v**2)) * 1 / sd1**2 - X1.T @ X1 / sd1**2
)
# add zeros for derv beta 0
der_b1_beta1 = np.concatenate(
(der_b1_beta1, np.zeros((n_col_X1, n_col_X0))), axis=1
)

# derivation wrt gamma
der_b1_gamma = -(X1Z1 * rho1v / (sd1 * (1 - rho1v ** 2)))
der_b1_gamma = -(X1Z1 * rho1v / (sd1 * (1 - rho1v**2)))
der_b1_gamma = np.concatenate((der_b1_beta1, der_b1_gamma), axis=1)
# derv wrt sigma 1
der_b1_sd = (
Expand All @@ -150,7 +151,7 @@ def calc_hess_beta1(
# derv wrt rho1
der_b1_rho = (
-(
eta1 * rho_aux1 * rho1v / ((1 - rho1v ** 2) ** 0.5)
eta1 * rho_aux1 * rho1v / ((1 - rho1v**2) ** 0.5)
+ norm.pdf(lambda1) / norm.cdf(lambda1) * rho_aux2
)
* 1
Expand All @@ -176,21 +177,21 @@ def calc_hess_beta0(
"""

# define some aux_vars
rho_aux1 = lambda0 * rho0v / (1 - rho0v ** 2) - nu0 / (1 - rho0v ** 2) ** 0.5
rho_aux2 = rho0v ** 2 / ((1 - rho0v ** 2) ** (3 / 2)) + 1 / (1 - rho0v ** 2) ** 0.5
sd_aux1 = rho0v ** 2 / (1 - rho0v ** 2)
sd_aux2 = rho0v / (np.sqrt(1 - rho0v ** 2))
rho_aux1 = lambda0 * rho0v / (1 - rho0v**2) - nu0 / (1 - rho0v**2) ** 0.5
rho_aux2 = rho0v**2 / ((1 - rho0v**2) ** (3 / 2)) + 1 / (1 - rho0v**2) ** 0.5
sd_aux1 = rho0v**2 / (1 - rho0v**2)
sd_aux2 = rho0v / (np.sqrt(1 - rho0v**2))

# add zeros for beta0
der_b0_beta1 = np.zeros((n_col_X1, n_col_X0))

# beta0
der_b0_beta0 = (
-(X0X0 * (rho0v ** 2 / (1 - rho0v ** 2)) * 1 / sd0 ** 2) + X0.T @ X0 / sd0 ** 2
-(X0X0 * (rho0v**2 / (1 - rho0v**2)) * 1 / sd0**2) + X0.T @ X0 / sd0**2
)
der_b0_beta0 = np.concatenate((der_b0_beta1, der_b0_beta0), axis=1)
# gamma
der_b0_gamma = -X0Z0 * rho0v / (1 - rho0v ** 2) * 1 / sd0
der_b0_gamma = -X0Z0 * rho0v / (1 - rho0v**2) * 1 / sd0
der_b0_gamma = np.concatenate((der_b0_beta0, der_b0_gamma), axis=1)

# add zeros for sigma1 and rho1
Expand All @@ -204,15 +205,15 @@ def calc_hess_beta0(
- 2 * nu0
)
* 1
/ sd0 ** 2
/ sd0**2
)
der_b0_sd = np.expand_dims((der_b0_sd.T @ X0), 1)
der_b0_sd = np.concatenate((der_b0_gamma, der_b0_sd), axis=1)

# rho
der_b0_rho = (
(
eta0 * -rho_aux1 * (rho0v / ((1 - rho0v ** 2) ** 0.5))
eta0 * -rho_aux1 * (rho0v / ((1 - rho0v**2) ** 0.5))
+ norm.pdf(lambda0) / (1 - norm.cdf(lambda0)) * rho_aux2
)
* 1
Expand Down Expand Up @@ -252,7 +253,7 @@ def calc_hess_gamma(

der_gamma_beta = np.zeros((Z1.shape[1], num_col_X1X0))

der_g_gamma = -(1 / (1 - rho1v ** 2) * Z1Z1 + 1 / (1 - rho0v ** 2) * Z0Z0)
der_g_gamma = -(1 / (1 - rho1v**2) * Z1Z1 + 1 / (1 - rho0v**2) * Z0Z0)
der_g_gamma = np.concatenate((der_gamma_beta, der_g_gamma), axis=1)

# sigma1
Expand All @@ -261,19 +262,19 @@ def calc_hess_gamma(
@ np.einsum("ij, i ->ij", X1, nu1)
/ sd1
* rho1v
/ (1 - rho1v ** 2)
/ (1 - rho1v**2)
)[:, 0]
der_g_sd1 = np.expand_dims(der_g_sd1, 0).T
der_g_sd1 = np.concatenate((der_g_gamma, der_g_sd1), axis=1)

# rho1
aux_rho11 = np.einsum("ij, i ->ij", Z1, eta1).T @ (
lambda1 * rho1v / (1 - rho1v ** 2) - nu1 / np.sqrt(1 - rho1v ** 2)
lambda1 * rho1v / (1 - rho1v**2) - nu1 / np.sqrt(1 - rho1v**2)
)
aux_rho21 = Z1.T @ (norm.pdf(lambda1) / norm.cdf(lambda1))

der_g_rho1 = -aux_rho11 * 1 / (np.sqrt(1 - rho1v ** 2)) - aux_rho21 * rho1v / (
(1 - rho1v ** 2) ** (3 / 2)
der_g_rho1 = -aux_rho11 * 1 / (np.sqrt(1 - rho1v**2)) - aux_rho21 * rho1v / (
(1 - rho1v**2) ** (3 / 2)
)
der_g_rho1 = np.expand_dims(der_g_rho1, 0).T
der_g_rho1 = np.concatenate((der_g_sd1, der_g_rho1), axis=1)
Expand All @@ -284,19 +285,19 @@ def calc_hess_gamma(
@ np.einsum("ij, i ->ij", X0, nu0)
/ sd0
* rho0v
/ (1 - rho0v ** 2)
/ (1 - rho0v**2)
)[:, 0]
der_g_sd0 = np.expand_dims(der_g_sd0, 0).T
der_g_sd0 = np.concatenate((der_g_rho1, -der_g_sd0), axis=1)

# rho1
aux_rho10 = np.einsum("ij, i ->ij", Z0, eta0).T @ (
lambda0 * rho0v / (1 - rho0v ** 2) - nu0 / np.sqrt(1 - rho0v ** 2)
lambda0 * rho0v / (1 - rho0v**2) - nu0 / np.sqrt(1 - rho0v**2)
)
aux_rho20 = -Z0.T @ (norm.pdf(lambda0) / (1 - norm.cdf(lambda0)))

der_g_rho0 = aux_rho10 * 1 / (np.sqrt(1 - rho0v ** 2)) + aux_rho20 * rho0v / (
(1 - rho0v ** 2) ** (3 / 2)
der_g_rho0 = aux_rho10 * 1 / (np.sqrt(1 - rho0v**2)) + aux_rho20 * rho0v / (
(1 - rho0v**2) ** (3 / 2)
)
der_g_rho0 = np.expand_dims(-der_g_rho0, 0).T
der_g_rho0 = np.concatenate((der_g_sd0, der_g_rho0), axis=1)
Expand Down Expand Up @@ -328,52 +329,52 @@ def calc_hess_dist(
Delta_sd1 = (
+1 / sd1
- (norm.pdf(lambda1) / norm.cdf(lambda1))
* (rho1v * nu1 / (np.sqrt(1 - rho1v ** 2) * sd1))
- nu1 ** 2 / sd1
* (rho1v * nu1 / (np.sqrt(1 - rho1v**2) * sd1))
- nu1**2 / sd1
)
Delta_sd1_der = (
nu1
/ sd1
* (
-eta1 * (rho1v ** 2 * nu1) / (1 - rho1v ** 2)
+ (norm.pdf(lambda1) / norm.cdf(lambda1)) * rho1v / np.sqrt(1 - rho1v ** 2)
-eta1 * (rho1v**2 * nu1) / (1 - rho1v**2)
+ (norm.pdf(lambda1) / norm.cdf(lambda1)) * rho1v / np.sqrt(1 - rho1v**2)
+ 2 * nu1
)
)

Delta_sd0 = (
+1 / sd0
+ (norm.pdf(lambda0) / (1 - norm.cdf(lambda0)))
* (rho0v * nu0 / (np.sqrt(1 - rho0v ** 2) * sd0))
- nu0 ** 2 / sd0
* (rho0v * nu0 / (np.sqrt(1 - rho0v**2) * sd0))
- nu0**2 / sd0
)
Delta_sd0_der = (
nu0
/ sd0
* (
-eta0 * (rho0v ** 2 * nu0) / (1 - rho0v ** 2)
-eta0 * (rho0v**2 * nu0) / (1 - rho0v**2)
- (norm.pdf(lambda0) / (1 - norm.cdf(lambda0)))
* rho0v
/ np.sqrt(1 - rho0v ** 2)
/ np.sqrt(1 - rho0v**2)
+ 2 * nu0
)
)

aux_rho11 = lambda1 * rho1v / (1 - rho1v ** 2) - nu1 / np.sqrt(1 - rho1v ** 2)
aux_rho12 = 1 / (1 - rho1v ** 2) ** (3 / 2)
aux_rho11 = lambda1 * rho1v / (1 - rho1v**2) - nu1 / np.sqrt(1 - rho1v**2)
aux_rho12 = 1 / (1 - rho1v**2) ** (3 / 2)

aux_rho_rho11 = (np.dot(gamma, Z1.T) * rho1v - nu1) / (1 - rho1v ** 2) ** (3 / 2)
aux_rho_rho11 = (np.dot(gamma, Z1.T) * rho1v - nu1) / (1 - rho1v**2) ** (3 / 2)
aux_rho_rho12 = (
2 * np.dot(gamma, Z1.T) * rho1v ** 2 + np.dot(gamma, Z1.T) - 3 * nu1 * rho1v
) / (1 - rho1v ** 2) ** (5 / 2)
2 * np.dot(gamma, Z1.T) * rho1v**2 + np.dot(gamma, Z1.T) - 3 * nu1 * rho1v
) / (1 - rho1v**2) ** (5 / 2)

aux_rho01 = lambda0 * rho0v / (1 - rho0v ** 2) - nu0 / np.sqrt(1 - rho0v ** 2)
aux_rho02 = 1 / (1 - rho0v ** 2) ** (3 / 2)
aux_rho01 = lambda0 * rho0v / (1 - rho0v**2) - nu0 / np.sqrt(1 - rho0v**2)
aux_rho02 = 1 / (1 - rho0v**2) ** (3 / 2)

aux_rho_rho01 = (np.dot(gamma, Z0.T) * rho0v - nu0) / (1 - rho0v ** 2) ** (3 / 2)
aux_rho_rho01 = (np.dot(gamma, Z0.T) * rho0v - nu0) / (1 - rho0v**2) ** (3 / 2)
aux_rho_rho02 = (
2 * np.dot(gamma, Z0.T) * rho0v ** 2 + np.dot(gamma, Z0.T) - 3 * nu0 * rho0v
) / (1 - rho0v ** 2) ** (5 / 2)
2 * np.dot(gamma, Z0.T) * rho0v**2 + np.dot(gamma, Z0.T) - 3 * nu0 * rho0v
) / (1 - rho0v**2) ** (5 / 2)

# for sigma1

Expand All @@ -385,7 +386,7 @@ def calc_hess_dist(
1
/ sd1
* (
-eta1 * aux_rho11 * (rho1v * nu1) / (np.sqrt(1 - rho1v ** 2))
-eta1 * aux_rho11 * (rho1v * nu1) / (np.sqrt(1 - rho1v**2))
- (norm.pdf(lambda1) / norm.cdf(lambda1)) * aux_rho12 * nu1
)
)
Expand All @@ -411,7 +412,7 @@ def calc_hess_dist(
1
/ sd0
* (
-eta0 * aux_rho01 * (rho0v * nu0) / (np.sqrt(1 - rho0v ** 2))
-eta0 * aux_rho01 * (rho0v * nu0) / (np.sqrt(1 - rho0v**2))
+ (norm.pdf(lambda0) / (1 - norm.cdf(lambda0))) * aux_rho02 * nu0
)
)
Expand Down
1 change: 1 addition & 0 deletions development/tests/property/property_auxiliary.py
@@ -1,4 +1,5 @@
"""This module contains some auxiliary functions for the property testing."""

import glob
import os

Expand Down
1 change: 1 addition & 0 deletions development/tests/property/run.py
@@ -1,4 +1,5 @@
"""The test provides the basic capabilities to run numerous property tests."""

import os

import functools
Expand Down
1 change: 1 addition & 0 deletions development/tests/reliability/reliability.py
Expand Up @@ -3,6 +3,7 @@
used. Additionally the module creates two different figures for the reliability section of the
documentation.
"""

import matplotlib.pyplot as plt
import numpy as np
import pandas as pd
Expand Down
1 change: 1 addition & 0 deletions development/tests/replication/run.py
@@ -1,6 +1,7 @@
"""This script replicates the estimation results from Cainero 2011 via the grmpy estimation method.
Additionally it returns a figure of the Marginal treatment effect based on the estimation results.
"""

import json
import matplotlib.pyplot as plt
import numpy as np
Expand Down
1 change: 1 addition & 0 deletions docs/__init__.py
@@ -1,4 +1,5 @@
"""The module allows to run tests from inside the interpreter."""

import os

import pytest
Expand Down
1 change: 1 addition & 0 deletions docs/source/figures/create.py
@@ -1,6 +1,7 @@
"""This module creates all available figures and then provides them in a compiled
document for review.
"""

import os
import subprocess

Expand Down
1 change: 1 addition & 0 deletions docs/source/figures/scripts/__init__.py
@@ -1,4 +1,5 @@
"""The module allows to run tests from inside the interpreter."""

import os

import pytest
Expand Down
1 change: 1 addition & 0 deletions docs/source/figures/scripts/fig-distribution-joint.py
@@ -1,4 +1,5 @@
"""This module creates the contour plots for the bivariate normal distribution."""

import numpy as np
import matplotlib.pyplot as plt
from scipy.stats import multivariate_normal
Expand Down
1 change: 1 addition & 0 deletions docs/source/figures/scripts/fig-eh-marginal-effect.py
@@ -1,6 +1,7 @@
"""The script creates a figure to illustrate the appearance of the marginal treatment
effect in the abscence and presence of individual heterogeneity.
"""

import numpy as np
import matplotlib.pyplot as plt

Expand Down
1 change: 1 addition & 0 deletions docs/source/figures/scripts/fig-local-average-treatment.py
@@ -1,5 +1,6 @@
"""This module contains the code for a local average treatment graph.
"""

import matplotlib.pyplot as plt

from fig_config import OUTPUT_DIR, RESOURCE_DIR
Expand Down
1 change: 1 addition & 0 deletions docs/source/figures/scripts/fig-reliability-par.py
@@ -1,6 +1,7 @@
"""This script replicates the results of Caneiro 2011 via using a
mock data set and plots the original as well as the estimated mar-
ginal treatment effect"""

import json

import pandas as pd
Expand Down
1 change: 1 addition & 0 deletions docs/source/figures/scripts/fig-treatment-effects.py
@@ -1,4 +1,5 @@
""""""

import numpy as np
import matplotlib.pyplot as plt
from scipy.stats import norm
Expand Down
1 change: 1 addition & 0 deletions docs/source/figures/scripts/fig-weights-marginal-effect.py
@@ -1,6 +1,7 @@
""" This script creates a figure to illustrate how the usual treatment effects can be
constructed by using differen weights on the marginal treatment effect.
"""

import numpy as np
import matplotlib.pyplot as plt
from scipy.stats import norm
Expand Down
1 change: 1 addition & 0 deletions docs/source/figures/scripts/fig_config.py
@@ -1,4 +1,5 @@
"""This file contains paths that are used for the execution of all figure scripts."""

import os

RESOURCE_DIR = os.path.dirname(os.path.abspath(__file__)) + "/resources"
Expand Down
1 change: 1 addition & 0 deletions grmpy/__init__.py
Expand Up @@ -8,6 +8,7 @@
gp.<func>
"""

import pytest

from grmpy.estimate.estimate import fit # noqa: F401
Expand Down
1 change: 1 addition & 0 deletions grmpy/estimate/estimate.py
@@ -1,6 +1,7 @@
"""The module provides an estimation process given the simulated data set and the
initialization file.
"""

from grmpy.check.auxiliary import read_data
from grmpy.check.check import (
check_est_init_dict,
Expand Down
1 change: 1 addition & 0 deletions grmpy/estimate/estimate_output.py
@@ -1,4 +1,5 @@
"""This module contains methods for producing the estimation output files."""

import time
from textwrap import wrap

Expand Down

0 comments on commit e89265f

Please sign in to comment.