text:
# -*- coding: utf-8 -*-
"""
Module used for performing Approximate Bayesian Computation (ABC).
Based on description in the Handbook of Approxiamte Bayesian Computation (2019).
"""

import os
from datetime import datetime
import warnings

import numpy as np
from scipy.optimize import minimize_scalar
from numba import jit

# custom modules: logging and storage
from . import log_simulations
from . import store_results

#custom modules: functions used in ABC
from . import distances
from . import kernels

########################### Function definitions ##############################
@jit(nopython = True)
def ESS(w):
    """
    Computes effective sample size based on a vector of weights
    (see: Handbook of Approxiamte Bayesian Computation (2019), p. 114).

    Arguments
    ---------
    w : ndarray
        1d vector of weights

    Returns
    --------
    ESS : float
        effective sample size
    """
    W = w/np.sum(w) #normalised weights
    return np.sum(W**2)**(-1)

def __default_timestamp_format__(datetime_object):
    return str(datetime_object).replace(":","-").split(".")[0]

########################### Class definitions #################################
class SMC_MCMC():
    """
    Class implementing the sequential markov-chain monte-carlo ABC algorithm.
    ========================================================================

    Following attributes are set during creation of a new instance.
    More attributes are set during run call.

    Attributes
    ----------
    generative_model : function
        generative model to be used for generating simulations.
    priors : function or container,
        Prior distributions from which the initial particle distribution is drawn.
        This can bne specified as:
            - a function that can draw samples from priors of the n_param parameters.
            - or a container (list, tuple) that unpacks to two ndarrays
            param_dist_0 and weights_0 of length n_param.
    n_param : int
        number of parameters to be inferred.
    y_obs : ndarray
        observed states (if s_obs == True, then y_obs is the observed summary statistics).
    summary_stat_f : function
        function that maps states y to summary statistics s.
    s_obs : bool, optional
        flag to determine if y_obs is s_obs instead. The default is False.
    n_particles : int, optional
        number of particles in the SMC distribution. The default is 100.
    n_draws_per_param : int, optional
        number of repeated simulations per parameter vector. The default is 1.
    fixed_model_parameters : dict, optional
        if set, these are the fixed arguments of generative_model. The default is None.
    zero_tol : float, optional
        zero tolerance. The default is 1e-10.
    max_SMC_iter : int, optional
        maximum number of SMC iterations. The default is 50.
    alpha_ESS : float, optional
        parameter to control the ESS decrease during every iteration.
        Value between 0 and 1.  The default is 0.9.
    ESS_min_factor : float, optional
        Minium accepted (relative) ESS before resampling . Value between 0 and 1. The default is 0.5.
    min_accept_ratio : float, optional
        Minimum ratio of moved particles vs candidates below which the algorithm is terminated.
        Value between 0 and 1. The default is 0.1.
    initial_bandwith : float, optional
        If set, this is the initial bandwith of the smoothing kernel.
        The default is None, in that case it is determined as 3 times
        the maximum (absolute) inial distance from s_obs.
    smoothing_kernel_f : str or function, optional
        smoothing kernel function to be used from the module kernels. The default is "uniform".
        Can also be a self-implemented function.
    distance_f : str or function, optional
        distance function to be used from the module kernels. The default is "Euclidean".
        Can also be a self-implemented function.
    MCMC_kernel_specifics : dict, optional
        if given, it's a dict containing the arguments to be passed the
        transition MCMC kernel object. The default is None.

    """
    def __init__(self, generative_model, priors, n_param,
                 y_obs,
                 summary_stat_f,
                 s_obs = False,
                 n_particles = 100,
                 n_draws_per_param = 1,
                 fixed_model_parameters = None,
                 zero_tol = 1e-10,
                 max_SMC_iter = 50,
                 alpha_ESS = 0.9,
                 ESS_min_factor = 0.5,
                 min_accept_ratio = 0,
                 initial_bandwith = None,
                 smoothing_kernel_f = "uniform",
                 distance_f = "Euclidean",
                 MCMC_kernel_specifics = None):

        # observations, model definition, priors
        #---------------------------------------
        self.summary_stat_f = summary_stat_f
        if s_obs:
            self.s_obs = y_obs
        else:
            self.s_obs = summary_stat_f(y_obs)

        self.dim_s = self.s_obs.shape
        if not self.dim_s: #empty tuple -> dimension 0
            self.dim_s = (1,)
        self.generative_model = generative_model
        if fixed_model_parameters is not None:
            self.fixed_model_parameters = fixed_model_parameters
        else:
            self.fixed_model_parameters = {}
        self.n_param = n_param
        self.priors = priors

        # hyperparameters
        #----------------
        self.initial_bandwith = initial_bandwith
        self.n_particles = n_particles
        self.n_draws_per_param = n_draws_per_param
        self.alpha_ESS = alpha_ESS
        self.zero_tol = zero_tol
        self.max_SMC_iter = max_SMC_iter
        self.ESS_min_factor = ESS_min_factor
        self.min_accept_ratio = min_accept_ratio

        # set functions to be used
        #-------------------------
        self.set_smoothing_kernel_f(smoothing_kernel_f)
        self.set_distance_f(distance_f)
        try:
            self.MCMC_kernel = MCMC_kernel(self, **MCMC_kernel_specifics)
        except TypeError:
            self.MCMC_kernel = MCMC_kernel(self) #use default

    def run(self, path_results = None, log_file = None, seed = None, **hyperparameters):
        """
        Run the SMC-MCMC algorithm for ABC inference.

        Parameters
        ----------
        path_results : str, optional
            Path to results folder. If None is passed, the results folder is
            `results_<timestamp>`. The default is None.
        log_file : str, optional
            Path to log file. If None is passed, the log file is
            ``. The default is None.
        seed : int, optional
            seed to be fed to the random number generator (`np.random.seed()`).
            The default is None.
        **hyperparameters : optional
            following hyperparameters can be set as keywords of run():
                initial_bandwith, n_particles, n_draws_per_param, alpha_ESS,
                zero_tol, max_SMC_iter, ESS_min_factor, min_accept_ratio, priors.
                For an explanation of these parameters: see help(SMC_MCMC).

        Returns
        -------
        param_distributions: ndarray,
            distributions of parameter vectors (particles) at every SMC iteration.

        weight_distributions: ndarray,
            particles weights at every SMC iteration.

        ABC_bandwiths: ndarray,
            smoothing kernel bandwith at every SMC iteration.

        """
        self.i_SMC = 0 # set counter of SMC iterations

        # (re)set hyperparameters
        #-----------------------------------------
        for hp_name, hp_val in hyperparameters.items():
            if hp_name in {"initial_bandwith",
                           "n_particles",
                           "n_draws_per_param",
                           "alpha_ESS",
                           "zero_tol",
                           "max_SMC_iter",
                           "ESS_min_factor",
                           "min_accept_ratio",
                           "priors"}:
                setattr(self,hp_name,hp_val)
            else:
                warnings.warn(f"Cannot set {hp_name} during run.")
        self.ESS_min = self.ESS_min_factor * self.n_particles

        # set time and seed
        #-----------------------------------------
        self.t_start = datetime.now()

        if seed is None:
            self.seed = int(self.t_start.timestamp()) # default: seed == timestamp
        else:
            self.seed = seed

        np.random.seed(self.seed)

        # initialise directories, files and arrays
        #-----------------------------------------
        try:
            store_results.mkdir_p(path_results)
            self.path_results = path_results
        except TypeError:
            self.path_results = f"results_{__default_timestamp_format__(self.t_start)}"
            store_results.mkdir_p(self.path_results)

        self.initialise_logger(log_file) #logger

        self.weight_distributions = np.empty((self.max_SMC_iter,self.n_particles),dtype = np.float64)
        self.param_distributions = np.empty((self.max_SMC_iter,self.n_particles, self.n_param),dtype = np.float64)
        self.ABC_bandwiths = np.empty(self.max_SMC_iter) #smoothing kernel bandwith at every iteration

        self.moved_particles = np.zeros(self.max_SMC_iter) # array for storing amount of particels that are moved at every SMC iteration
        self.move_candidates = self.moved_particles.copy() # array for storing amount of particles that were selected to move at every SMC iteration

        # assign starting values
        #-----------------------
        try:
            self.param_dist_0, self.weights_0 = self.draw_initial_param_dist(self.n_particles, self.n_param, self.priors)
        except TypeError:
            self.param_dist_0, self.weights_0 = self.priors
        self.param_distributions[0] = self.param_dist_0 # add to results array
        self.weight_distributions[0] = self.weights_0   # add to results array
        self.MCMC_kernel.update_scale("dynamic_param_dist_range") # update kernel if dynamic_param_dist_range mode is used

        self.d_0 = self.generate_distances_to_obs(self.param_dist_0)
        if self.initial_bandwith is None:
            self.initial_bandwith = np.max(np.abs(self.d_0))*3 # initial smoothing kernel bandwith
        self.ABC_bandwiths[0] = self.initial_bandwith # add to results array (-> see: initialise_results)

        self.logger_header.info("iter\t bandwith\tESS")
        self.logger_header.info(f"{0:3d}\t {self.initial_bandwith:>10.2f}\t{ESS(self.weight_distributions[0]):>3.0f}\t")

        # Start iteration
        #======================================================================
        for i_SMC in range(1,self.max_SMC_iter):
            self.i_SMC = i_SMC # set counter of SMC iterations
            self.log_message = "" # will be extended during iteration

            # step 1: reweight
            #-----------------------------------------------------------------
            warnings.filterwarnings('error', category= RuntimeWarning) # catch runtimewarnings
            try:
                self.ABC_bandwiths[i_SMC], weights_1 = self.reweight(self.weights_0, self.ABC_bandwiths[i_SMC-1], self.d_0 )
                self.log_message += f"{i_SMC:3d}\t {self.ABC_bandwiths[i_SMC]:>10.2f}\t{ESS(weights_1):>3.0f}: "
            except RuntimeWarning as rn:
                # break if (one or more) weights blow up
                self.logger_header.exception(f"""
-------------------------------------------------------------------------------------------------------------------------
Iteration {i_SMC:3d}:
     runtimewarning {rn} during reweight.
                                          """,exc_info = True)
                self.logger_header.info(f"""SMC-MCMC algorithm stopped at                  {datetime.now()}
=========================================================================================================================""")
                warnings.filterwarnings('once', category= RuntimeWarning)
                return self.param_distributions[:i_SMC], self.weight_distributions[:i_SMC], self.ABC_bandwiths[:i_SMC]
            except AssertionError:
                self.logger_header.exception(f"""
-------------------------------------------------------------------------------------------------------------------------
Iteration {i_SMC:3d}:
    zero in kernel ratio denominator during reweight.
                                          """,exc_info = True)
                self.logger_header.info(f"""SMC-MCMC algorithm stopped at                  {datetime.now()}
=========================================================================================================================""")
                warnings.filterwarnings('once', category= RuntimeWarning)
                return self.param_distributions[:i_SMC], self.weight_distributions[:i_SMC], self.ABC_bandwiths[:i_SMC]
            warnings.filterwarnings('once', category= RuntimeWarning)

            if np.sum(weights_1) == 0: # check if all weights are more than 0
                self.logger_header.info(self.log_message)
                self.logger_header.exception(f"""
-------------------------------------------------------------------------------------------------------------------------
Iteration {i_SMC:3d}:
    all particle weights are zero.

SMC-MCMC algorithm stopped at                  {datetime.now()}
=========================================================================================================================""")
                return self.param_distributions[:i_SMC], self.weight_distributions[:i_SMC], self.ABC_bandwiths[:i_SMC]

            # step 2: resample (if necessary)
            #-----------------------------------------------------------------
            weights_1_norm = weights_1/np.sum(weights_1) # normalize weights
            if ESS(weights_1) < self.ESS_min:
                self.log_message+="(resample) "
                weights_1, weights_1_norm = self.resample(weights_1, weights_1_norm)
                self.MCMC_kernel.update_scale("dynamic_param_dist_range") # update kernel if dynamic_param_dist_range mode is used
            else:
                self.log_message+="           "

            self.weights_1 = weights_1

            # step 3: move particles
            #-----------------------------------------------------------------
            self.candidate_i, self.param_perturbed, self.d_perturbed = self.perturb(self.param_dist_0, weights_1_norm )
            n_candidates = np.sum(self.candidate_i)
            warnings.filterwarnings('error', category= RuntimeWarning)
            try:
                kernel_ratio = self.calculate_kernel_ratio(d_0 = self.d_0[self.candidate_i],
                                                           d_1 = self.d_perturbed,
                                                           eps_0 = self.ABC_bandwiths[i_SMC],
                                                           eps_1 = self.ABC_bandwiths[i_SMC])
            except RuntimeWarning as rn:
                self.logger_header.exception(f"""
-------------------------------------------------------------------------------------------------------------------------
Iteration {i_SMC:3d}:
     runtimewarning {rn} during moving.
                                          """,exc_info = True)
                self.logger_header.info(f"""SMC-MCMC algorithm stopped at                  {datetime.now()}
=========================================================================================================================""")
                warnings.filterwarnings('once', category= RuntimeWarning)
                return self.param_distributions[:i_SMC], self.weight_distributions[:i_SMC], self.ABC_bandwiths[:i_SMC]
            except AssertionError:
                self.logger_header.exception(f"""
-------------------------------------------------------------------------------------------------------------------------
Iteration {i_SMC:3d}:
    zero in kernel ratio denominator during moving.
                                          """,exc_info = True)
                self.logger_header.info(f"""SMC-MCMC algorithm stopped at                  {datetime.now()}
=========================================================================================================================""")
                warnings.filterwarnings('once', category= RuntimeWarning)
                return self.param_distributions[:i_SMC], self.weight_distributions[:i_SMC], self.ABC_bandwiths[:i_SMC]
            warnings.filterwarnings('once', category= RuntimeWarning)
            self.moved_i, self.moved_perturbed_i = self.accept_perturbed_particles(kernel_ratio, self.candidate_i)
            n_moved = np.sum(self.moved_perturbed_i)

            # determine new param dist

            self.param_dist_0[self.moved_i] = self.param_perturbed[self.moved_perturbed_i]
            self.weights_0 = self.weights_1
            self.d_0[self.moved_i] = self.d_perturbed[self.moved_perturbed_i]

            # add to results
            self.param_distributions[i_SMC] = self.param_dist_0.copy()
            self.weight_distributions[i_SMC] = self.weights_0.copy()
            self.move_candidates[i_SMC] = n_candidates
            self.moved_particles[i_SMC] = n_moved

            self.log_message += f"{n_moved:4d} particles moved of the {n_candidates:4d} candidates"

            if n_moved/n_candidates < self.min_accept_ratio:
                self.log_message += f"""
-------------------------------------------------------------------------------------------------------------------------
Iteration {i_SMC:3d}: accept ratio dropped below minimum value

SMC-MCMC algorithm stopped at                  {datetime.now()}
========================================================================================================================="""

                self.logger_header.info(self.log_message)
                return self.param_distributions[:i_SMC+1], self.weight_distributions[:i_SMC+1], self.ABC_bandwiths[:i_SMC+1]
            else:
                self.logger_header.info(self.log_message)

        # stop when maximum of iterations is reached
        self.logger_header.info(f"""
-------------------------------------------------------------------------------------------------------------------------
Max number of {self.max_SMC_iter-1} iterations reached

SMC-MCMC algorithm stopped at                  {datetime.now()}
=========================================================================================================================""")
        return self.param_distributions, self.weight_distributions, self.ABC_bandwiths

    def set_smoothing_kernel_f(self, kernel_f):
        """
        Sets the smoothing kernel function as of SMC_MCMC object.

        Parameters
        ----------
        kernel_f : str or function
            kernel function to be set. This can be one defined in the kernels
            module or can be self-implemented.

        Sets attribute
        --------------
        smoothing_kernel_f : function
            smoothing kernel function (used in ABC approximation)

        """
        try:
            self.smoothing_kernel_f = getattr(kernels,kernel_f)
        except TypeError:
            self.smoothing_kernel_f = kernel_f

    def set_distance_f(self, distance_f):
        """
        Sets the distance function as method of SMC_MCMC object.

        Parameters
        ----------
        distance_f :
            distance function to be set.  This can be one defined in the distances
            module or can be self-implemented.

        """
        try:
            self.distance_f = getattr(distances,distance_f)
        except TypeError:
            self.distance_f = distance_f

    def initialise_logger(self, log_file):
        """
        Initialise logging.

        Two logger objects (module: logging) are created and set as attributes.

        Parameters
        ----------
        log_file : str, optional
            path to log file. If None is passed, the log file is
            `<self.path_results>/log_SMC_ABC_<timestamp>.log`. The default is None.

        Returns
        -------
        logger_header : logger object
            logger used as header in the log file

        """

        if log_file is None:
            timestamp_formatted = __default_timestamp_format__(self.t_start)
            log_file = f"{self.path_results}/log_SMC_ABC_{timestamp_formatted}.log"
        branch, commit_name, commit_hash = log_simulations.get_last_git_commit()

        # create all directories in path to log file
        store_results.mkdir_p(os.path.split(log_file)[0])

        self.logger_header = log_simulations.initialise_logger("", log_file, file_format_str = '%(message)s')

        # initialise log file with following header
        self.logger_header.info(f"""
=========================================================================================================================
Approximate Bayesian Computation

Version:                                       [{branch}]-{commit_name}-{commit_hash}
SMC-MCMC algorithm started at                  {self.t_start}
with seed {self.seed}
=========================================================================================================================
                                """)

    def generate_distances_to_obs(self, param_dist):
        """
        Generate simulated summary statistics from `param_dist` and claculate
        distance to `s_obs`.

        Parameters
        ----------
        param_dist : ndarray
            Parameter distribution , of which all the particles (parameter vectors)
            are fed to the generative model.
            shape: (n_particles,n_param)

        Returns
        -------
        distances: ndarray
            Distances between simulated and observed summary statistics.
            shape: (n_particles,)

        """
        s_param_dist = self.generative_model(param_dist,self.n_draws_per_param,
                                             self.summary_stat_f,
                                             self.dim_s,
                                             **self.fixed_model_parameters)
        return self.distance_f(s_param_dist,self.s_obs)

    def draw_initial_param_dist(self, n_particles, n_param, priors):
        """
        Draws an initial parameter distribution from priors and sets this as
        attribute of SMC_MCMC object.

        Returns
        ---------------
        param_dist_0 : ndarray
            Initial parameter distribution, drawn from prior.
            shape: (n_particles,n_param)
        weights_0 : ndarray
            Initial weights distribution (all equal to 1/n_particles)
        """
        param_dist_0 = priors(n_particles).reshape(n_particles,n_param) # draws from priors
        weights_0 = np.full((n_particles), 1/n_particles )# initial weights
        return param_dist_0, weights_0

    def reweight(self, weights, eps_0, distances):
        """
        Reweight step:
            - find kernel kernel bandwith at which ESS(w_0, eps_1) == alpha*ESS(w_0, eps_0)
            - update weights based on these new weights

        Returns
        -------
        w_1 : ndarray
            weights obtained after adjusting the smoothing kernel bandwith

        Returns
        -------
        eps_1: float,
            new smoothing kernel bandwith.

        w_1: ndarray,
            new particle weights.

        """
        opt_bandwith = minimize_scalar(self.ESS_objective_fun,
                                        args = (eps_0, weights, distances),
                                        method='bounded',
                                        bounds = (0, eps_0))
        eps_1 = opt_bandwith.x

        # compute new weights
        w_1 = self.update_weights(w_0 = weights,
                                  d_0 = distances,
                                  d_1 = distances,
                                  eps_0 = eps_0,
                                  eps_1 = eps_1)
        return eps_1, w_1

    def resample(self, w_1, w_1_norm):
        """
        Perform resample step.

        Parameters
        ----------
        w_1 : ndarray
            particle weights.
        w_1_norm : ndarray
            normalised weights.

        Returns
        -------
        w_1 : ndarray
            particle weights.
        w_1_norm : ndarray
            normalised weights.

        Updates attributes
        ------------------
        param_dist_0 : ndarray,
            current particle distribution

        d_0 : ndarray,
            distances of simulations with current particles to s_obs

        """

        # select particle indices with probability == normalised weight
        param_dist_0_i = np.random.choice(np.arange(self.n_particles),
                                          size = self.n_particles,
                                          p = w_1_norm)
        self.param_dist_0 = self.param_dist_0[param_dist_0_i,]
        w_1 = w_1_norm = np.full((self.n_particles), 1/self.n_particles) # set all particle weights to be equal
        self.d_0 = self.d_0[param_dist_0_i,]
        return w_1, w_1_norm

    def perturb(self, param_dist_0, w_1_norm):
        """
        Perturb parameter vectors of a given parameter distribution.

        Parameters
        ----------
        param_dist_0 : ndarray
            current particle distribution.
        w_1_norm : ndarray
            normalised weights.

        Returns
        -------
        candidate_i : ndarray
            boolean indexes of candidates for move step.
        param_perturbed : ndarray
            perturbed parameters at positions specified in candidate_i param_dist_0[candidate_i].
        d_perturbed : ndarray
            distances of simulations with perturbed particles to s_obs (only updated if resampling was performed).

        """
        # select candidate particles for moving (normalised weights > w_tol)
        candidate_i = ~np.isclose(w_1_norm, 0, atol = self.zero_tol)

        # move candidates and calculate distances
        param_perturbed = self.MCMC_kernel.perturb_particles(param_dist_0[candidate_i,])
        d_perturbed = self.generate_distances_to_obs(param_perturbed)

        return candidate_i, param_perturbed, d_perturbed

    def accept_perturbed_particles(self, kernel_ratio, candidate_i):
        p_move = np.min(
                    np.pad(kernel_ratio.reshape(-1,1),
                           ((0,0),(1,0)),mode="constant", constant_values = 1),
                           axis = 1) # Handbook of Approxiamte Bayesian Computation (2019), p. 115

        ## pick moved candidates based on acceptance probability
        moved_perturbed_i = p_move > np.random.sample(size = p_move.shape)
        moved_i = np.zeros(candidate_i.shape,dtype=bool)
        moved_i[candidate_i] = moved_perturbed_i
        return moved_i, moved_perturbed_i

    def update_weights(self, w_0, d_0, d_1, eps_0, eps_1):
        """
        Update particle weights based on MCMC.

        Weights are based on the distances between the `n_draws_per_param` replicates
        of simulations, obtained  by feeding `n_particles` samples (particles)
        of the empirical parameter distribution to the generative model.

        Parameters
        ----------
        w_0 : ndarray
            1d vector of original weights
        d_0, d_1 : ndarray
            Arrays of (n_particles, n_draws_per_param) with original and new distances.
        eps_0, eps_1 : float
            kernel bandwiths

        Returns
        -------
        w_1 : ndarray
            updated weights.
        """

        w_1 = np.zeros_like(w_0)

        # exclude zero weights
        i_nonzero = np.logical_not(np.isclose(w_0, 0))
        w_1[i_nonzero] = w_0[i_nonzero]*self.calculate_kernel_ratio(d_0[i_nonzero], d_1[i_nonzero], eps_0, eps_1) # reweight
        return w_1

    def calculate_kernel_ratio(self,
                               d_0, d_1,
                               eps_0, eps_1):
        """
        Computes the kernel ratio, used in the MCMC steps.
        (see: Handbook of Approxiamte Bayesian Computation (2019), p. 114-115).

        Parameters
        ----------
        d_0, d_1 : ndarray
            Arrays of shape (n_particles, n_draws_per_param) with original and new distances.

        eps_0, eps_1 : float
            kernel bandwiths

        Returns
        -------
        kernel_ratios : ndarray
            vector of kernel ratios (values below zero_tol are set to zero)

        """
        # check if denominator is not zero
        denom = np.sum(self.smoothing_kernel_f(d_0, eps_0),axis = 1)
        try:
            assert not np.any(denom == 0), "zero in denominator of kernel ratio"
        except AssertionError:
            np.save(os.path.join(self.path_results,"kernel_ratio_error_denom.npy"), denom)
            np.save(os.path.join(self.path_results,"kernel_ratio_error_d_0.npy"), denom)
            np.save(os.path.join(self.path_results,"kernel_ratio_error_eps_0.npy"), denom)
            raise
        # compute kernel ratios
        kernel_ratios = np.sum(self.smoothing_kernel_f(d_1,eps_1),axis = 1)/denom
        return kernel_ratios

    def ESS_objective_fun(self, eps_1, eps_0, w_0, d):
            """
            Objective function to find the bandwith of the smoothing kernel K_eps for
            which the effective sample size of the weights w_1 is alpha*w_0,
            where w_0 is the weight vector that was calculated using the bandwith eps_0
            and distance matrix d_0_t at iteration m-1. w_1 are reweighted based on  the original weights.

            Arguments
            ---------
            eps_1 : scalar (float),
                new kernel bandwith

            eps_0 : scalar (float),
                original kernel bandwith

            w_0 : ndarray,
                1d vector of original weights

            d : ndarray,
                (N_pop x T) distance matrix of T distances obtained for every out of N_pop particles

            Returns
            --------
            out: scalar,
                |ESS(w_1) - alpha*ESS(w_0)|
            """
            w_1 = self.update_weights(w_0, d, d, eps_0, eps_1)
            return abs(ESS(w_1) - self.alpha_ESS*ESS(w_0))

class MCMC_kernel():
    """
    Class implementing an MCMC kernel.
    ========================================================================

    Attributes
    ----------

    MCMC_object : MCMC_object
        ABC object in which the MCMC kernel is used

    scale : float
        scale parameter (or bandwith) of kernel

    scale_mode : str
        method in which the kernel scale is updated.
        - "constant" : scale stays constant during iterations
        - "dynamic_param_dist_range" : scale is set as a fixed ratio (`r_MCMC`) of parameter range at every resampling.

    kernel_f : str
         type of kernel function to be used
         Supported types are "gaussian" and "uniform"

    **kwargs :
        optional keyword arguments to be set as attributes.
    """
    def __init__(self, MCMC_object, scale = 0.1, scale_mode = "constant",
                 kernel_f = "gaussian", bounds = None, **kwargs):
        self.MCMC_object = MCMC_object
        self.scale = scale
        self.scale_mode = scale_mode
        if kernel_f not in {"gaussian", "uniform"}:
            raise ValueError(f"Unsupported kernel function {kernel_f}")
        self.kernel_f = kernel_f
        if bounds is None or bounds.shape == (2,MCMC_object.n_param):
            self.bounds = bounds
        for kw, val in kwargs.items():
            setattr(self,kw,val)

    def update_scale(self, mode):
        if mode == self.scale_mode == "dynamic_param_dist_range":
            # scale of MCMC kernel is proportional to current range of param dist
            self.scale = np.ptp(self.MCMC_object.param_dist_0, axis = 0) * self.r_MCMC

    def perturb_particles(self, param_dist):
        """

        Parameters
        ----------
        param_dist : ndarray
            particles to be perturbed by the MCMC kernel, within parameter bounds (if these are given). shape (n_perturbed, n_param).
.

        Returns
        -------
        perturbed_param_dist : ndarray
            particles perturbed by the MCMC kernel, within parameter bounds (if these are given). shape (n_perturbed, n_param).

        """
        if self.kernel_f == "gaussian":
            perturbation = np.random.normal(scale = self.scale/2, size = param_dist.shape)
        elif self.kernel_f == "uniform":
            perturbation = np.random.uniform(low = -self.scale/2, high = self.scale/2 , size = param_dist.shape)

        perturbed_param_dist = param_dist + perturbation

        if self.bounds is None:
            return perturbed_param_dist

        # keep parameter values within bounds
        low_bounds, high_bounds = self.bounds[0,:], self.bounds[1,:]
        self.check_overshoot(perturbation, low_bounds, high_bounds)
        return self.reflect_at_boundaries(perturbed_param_dist, low_bounds, high_bounds)

    def check_overshoot (self, perturbation, low_bounds, high_bounds):
        """
        Check if perturbation doesn't overshoot the parameter range.

        Parameters
        ----------
        perturbation : ndarray
            perturbation by the MCMC kernel. shape (n_perturbed, n_param)
        low_bounds : ndarray
            lower bounds of parameters (if unbounded: nan).
        high_bounds : ndarray
            higher bounds of parameters (if unbounded: nan).

        """
        param_range = high_bounds - low_bounds
        i_not_nan = ~np.isnan(param_range)
        assert not np.any(np.abs(perturbation[:,i_not_nan]) > param_range[i_not_nan] ), "MCMC kernel overshoots parameter range"

    def reflect_at_boundaries(self, perturbed_param_dist, low_bounds, high_bounds):
        """
        Keep parameters within bounds.

        Parameter particles are reflected when they cross given boundaries.

        Parameters
        ----------
        perturbed_param_dist : ndarray
            particles perturbed by the MCMC kernel. shape (n_perturbed, n_param)
        low_bounds : ndarray
            lower bounds of parameters (if unbounded: nan).
        high_bounds : ndarray
            higher bounds of parameters (if unbounded: nan).

        Returns
        -------
        perturbed_param_dist : ndarray
            particles perturbed by the MCMC kernel, within parameter bounds. shape (n_perturbed, n_param)

        """
        for i in range(self.MCMC_object.n_param):
            if not np.isnan(low_bounds[i]):
                too_low = perturbed_param_dist[:,i] < low_bounds[i]
                perturbed_param_dist[too_low,i] = 2*low_bounds[i] - perturbed_param_dist[too_low,i]
            if not np.isnan(high_bounds[i]):
                too_high = perturbed_param_dist[:,i] > high_bounds[i]
                perturbed_param_dist[too_high,i] = 2*high_bounds[i] - perturbed_param_dist[too_high,i]
        return perturbed_param_dist

reason: import numpy,from scipy,from numba


text:
import numpy as np
from copy import deepcopy
from models.multiple_solution.root_multiple import RootAlgo

class BaseBMO(RootAlgo):
    """
    The blue monkey: A new nature inspired metaheuristic optimization algorithm
    """
    ID_POS = 0      # position
    ID_FIT = 1      # fitness
    ID_RAT = 2      # rate
    ID_WEI = 3      # weight

    def __init__(self, root_algo_paras=None, bmo_paras=None):
        RootAlgo.__init__(self, root_algo_paras)
        self.epoch = bmo_paras["epoch"]
        self.pop_size = bmo_paras["pop_size"]
        self.bm_teams = bmo_paras["bm_teams"]       # Number of blue monkey teams (5, 10, 20, ...)
        self.bm_size = int(self.pop_size/2)         # Number of all blue monkey
        self.bm_numbers = int(self.bm_size / self.bm_teams)     # Number of blue monkey in each team

    def _create_solution__(self, minmax=0):
        solution = np.random.uniform(self.domain_range[0], self.domain_range[1], self.problem_size)
        fitness = self._fitness_model__(solution=solution, minmax=minmax)
        rate = np.random.uniform(0, 1)
        weight = np.random.uniform(4, 6)
        return [solution, fitness, rate, weight]

    def _create_population__(self):
        t1 = []
        for i in range(self.bm_size):
            t2 = [self._create_solution__(self.ID_MIN_PROBLEM) for _ in range(self.bm_numbers)]
            t1.append(t2)
        t2 = [self._create_solution__(self.ID_MIN_PROBLEM) for _ in range(self.bm_size)]
        return t1, t2

    def _train__(self):
        bm_pop, child_pop = self._create_population__()

        best = []
        for items in bm_pop:
            bt = self._get_global_best__(items, self.ID_FIT, self.ID_MIN_PROBLEM)
            best.append(deepcopy(bt))
        g_best = deepcopy(self._get_global_best__(best, self.ID_FIT, self.ID_MIN_PROBLEM))
        child_best = self._get_global_best__(child_pop, self.ID_FIT, self.ID_MIN_PROBLEM)
        if g_best[self.ID_FIT] > child_best[self.ID_FIT]:
            g_best = deepcopy(child_best)

        for epoch in range(self.epoch):

            for items in bm_pop:
                items = sorted(items, key=lambda temp: temp[self.ID_FIT])
                items[-1] = deepcopy(child_best)

            for items in bm_pop:
                leader = self._get_global_best__(items, self.ID_FIT, self.ID_MIN_PROBLEM)
                for i in range(self.bm_numbers):
                    rate = 0.7 * items[i][self.ID_RAT] + (leader[self.ID_WEI] - items[i][self.ID_WEI]) * \
                              np.random.uniform() * (leader[self.ID_POS] - items[i][self.ID_POS])
                    pos = items[i][self.ID_POS] + np.random.uniform() * rate
                    pos = self._amend_solution_and_return__(pos)
                    fit = self._fitness_model__(pos)
                    if fit < items[i][self.ID_FIT]:
                        we = items[i][self.ID_WEI]
                        items[i] = [pos, fit, rate, we]

            for i in range(self.bm_size):
                rate = 0.7 * child_pop[i][self.ID_RAT] + (child_best[self.ID_WEI] - child_pop[i][self.ID_WEI]) * \
                    np.random.uniform() * (child_best[self.ID_POS] - child_pop[i][self.ID_POS])
                pos = child_pop[i][self.ID_POS] + np.random.uniform() * rate
                pos = self._amend_solution_and_return__(pos)
                fit = self._fitness_model__(pos)
                if fit < child_pop[i][self.ID_FIT]:
                    we = child_pop[i][self.ID_WEI]
                    child_pop[i] = [pos, fit, rate, we]

            current_child_best = self._get_global_best__(child_pop, self.ID_FIT, self.ID_MIN_PROBLEM)
            if current_child_best[self.ID_FIT] < child_best[self.ID_FIT]:
                child_best = deepcopy(current_child_best)
                if current_child_best[self.ID_FIT] < g_best[self.ID_FIT]:
                    g_best = deepcopy(current_child_best)

            for i in range(self.bm_teams):
                bt = self._get_global_best__(bm_pop[i], self.ID_FIT, self.ID_MIN_PROBLEM)
                best[i] = bt
            current_g_best = deepcopy(self._get_global_best__(best, self.ID_FIT, self.ID_MIN_PROBLEM))
            if current_g_best[self.ID_FIT] < g_best[self.ID_FIT]:
                g_best = deepcopy(current_g_best)

            self.loss_train.append(g_best[self.ID_FIT])
            if self.print_train:
                print("Epoch = {}, Fit = {}".format(epoch + 1, g_best[self.ID_FIT]))
        return g_best[self.ID_FIT], self.loss_train

ext: py

max_issues_repo_name: luanedge/metaheuristics

max_forks_repo_name: luanedge/metaheuristics

reason: import numpy


text:
import numpy as np
import tensorflow as tf
import matplotlib.pyplot as plt
import cv2

cascade_path='/home/alex/Downloads/haarcascade_frontalface_default.xml'
model_paramaters='/home/alex/Downloads/face_mask1.h5'
name="face_mask_detection"
new_model=tf.keras.models.load_model(model_paramaters)
cv2.namedWindow(name)
cap=cv2.VideoCapture(0)

while True:
	image_size=224
	_,frame= cap.read()
	faceCascade=cv2.CascadeClassifier(cascade_path)
	gray=cv2.cvtColor(frame, cv2.COLOR_BGR2GRAY)
	faces=faceCascade.detectMultiScale(gray,1.3,5)
	for x,y,w, h in faces:
		#burada sadece insan yüzünün oldugu kisimin matrix alinir
		roi_face=frame[y:y+h,y:x+w]
		final_array=cv2.resize(roi_face,(image_size,image_size))
		final_array=np.expand_dims(final_array,axis=0)
		final_array=final_array/255
		predictions=new_model.predict(final_array)
		font=cv2.FONT_HERSHEY_PLAIN
		font_scale=1.5
		print("predictions",predictions)
		#Maskeli durum
		if predictions<0.5:
			status='Mask {} %'.format(round(100-predictions[0][0]*100,2))
			cv2.rectangle(frame,(x,y),(x+w,y+h),(0,255,0),3)
			cv2.putText(frame,status,(x,y-45), font,font_scale,color=(0,255,0),thickness=3)
		else:
			status='No Mask {} %'.format(round(predictions[0][0]*100,2))
			cv2.rectangle(frame,(x,y),(x+w,y+h),(0,0,255),3)
			cv2.putText(frame,status,(x,y-45), font,font_scale,color=(0,0,255),thickness=3)

	cv2.imshow('face_mask_detection',frame)

	if cv2.waitKey(1)==27:
		break
cap.release()
cv2.destroyWindow()


text:
# -*- coding: utf-8 -*-
"""
@author: clausmichele
"""

import time
import tensorflow as tf
import numpy as np
import os
from utilis import *
from tensorflow.python.tools.inspect_checkpoint import print_tensors_in_checkpoint_file

def Temp3CNN(input, is_training=True, output_channels=3, reuse=tf.AUTO_REUSE):
	input_middle = input[:,:,:,3:6]
	with tf.variable_scope('temp-block1',reuse=reuse):
		output = tf.layers.conv2d(input, 128, 3, padding='same', activation=tf.nn.leaky_relu)
	for layers in range(2, 20):
		with tf.variable_scope('temp-block%d' % layers,reuse=reuse):
			output = tf.layers.conv2d(output, 64, 3, padding='same', name='conv%d' % layers, use_bias=False)
			output = tf.nn.leaky_relu(output)
	with tf.variable_scope('temp-block20', reuse=reuse):
		output = tf.layers.conv2d(output, output_channels, 3, padding='same', use_bias=False)
	return input_middle - output

def shuffle_in_unison(a, b):
	rng_state = np.random.get_state()
	np.random.shuffle(a)
	np.random.set_state(rng_state)
	np.random.shuffle(b)

class TemporalDenoiser(object):
	def __init__(self, sess, input_c_dim=9, batch_size=32):
		self.sess = sess
		self.input_c_dim = input_c_dim
		# build model
		self.Y_ = tf.placeholder(tf.float32, [None, None, None, self.input_c_dim],name='clean_frames')
		self.is_training = tf.placeholder(tf.bool, name='is_training')
		self.X = tf.placeholder(tf.float32, [None, None, None, self.input_c_dim],name='noisy_frames')
		self.Y = Temp3CNN(self.X, is_training=self.is_training)
		self.loss = (1.0 / batch_size) * tf.nn.l2_loss(self.Y_[:,:,:,3:6] - self.Y)
		self.lr = tf.placeholder(tf.float32, name='learning_rate')
		self.eva_psnr = tf_psnr(self.Y, self.Y_[:,:,:,3:6])
		optimizer = tf.train.AdamOptimizer(self.lr, name='AdamOptimizer')
		update_ops = tf.get_collection(tf.GraphKeys.UPDATE_OPS)
		with tf.control_dependencies(update_ops):
			self.train_op = optimizer.minimize(self.loss)
		init = tf.global_variables_initializer()
		self.sess.run(init)
		print("[*] Initialize model successfully...")

	def train(self, data, data_noisy, eval_data, eval_data_noisy, batch_size, ckpt_dir, epoch, lr, eval_every_epoch=2):
		# assert data range is between 0 and 1
		numBatch = int(data.shape[0] / batch_size)
		# load pretrained model
		load_model_status, global_step = self.load(ckpt_dir)
		if load_model_status:
			iter_num = global_step
			start_epoch = global_step // numBatch
			start_step = global_step % numBatch
			print("[*] Model restore successfully!")
		else:
			iter_num = 0
			start_epoch = 0
			start_step = 0
			print("[*] No pretrained model found!")
		# make summary
		tf.summary.scalar('loss', self.loss)
		tf.summary.scalar('lr', self.lr)
		writer = tf.summary.FileWriter('./logs', self.sess.graph)
		merged = tf.summary.merge_all()
		print("[*] Start training, with start epoch %d start iter %d : " % (start_epoch, iter_num))
		start_time = time.time()
		for epoch in range(start_epoch, epoch):
			shuffle_in_unison(data,data_noisy)
			for batch_id in range(start_step, numBatch):
				batch_images = data[batch_id * batch_size:(batch_id + 1) * batch_size, :, :, :]
				batch_images = batch_images.astype(np.float32) / 255.0 # normalize the data to 0-1
				batch_noisy = data_noisy[batch_id * batch_size:(batch_id + 1) * batch_size, :, :, :]
				batch_noisy = batch_noisy.astype(np.float32) / 255.0 # normalize the data to 0-1

				_, loss, summary = self.sess.run([self.train_op, self.loss, merged],
												 feed_dict={self.Y_: batch_images, self.X: batch_noisy, self.lr: lr[epoch],
															self.is_training: True})
				print("Epoch: [%2d] [%4d/%4d] time: %4.4f, loss: %.6f"
					  % (epoch + 1, batch_id + 1, numBatch, time.time() - start_time, loss))
				iter_num += 1
				writer.add_summary(summary, iter_num)
			self.save(iter_num, ckpt_dir)
		print("[*] Finish training.")

	def save(self, iter_num, ckpt_dir, model_name='Temp3-CNN'):
		saver = tf.train.Saver()
		checkpoint_dir = ckpt_dir
		if not os.path.exists(checkpoint_dir):
			os.makedirs(checkpoint_dir)
		print("[*] Saving model...")
		saver.save(self.sess,os.path.join(checkpoint_dir, model_name),global_step=iter_num)

	def load(self, checkpoint_dir):
		print("[*] Reading checkpoint...")
		saver = tf.train.Saver()
		ckpt = tf.train.get_checkpoint_state(checkpoint_dir)
		if ckpt and ckpt.model_checkpoint_path:
			full_path = tf.train.latest_checkpoint(checkpoint_dir)
			global_step = int(full_path.split('/')[-1].split('-')[-1])
			saver.restore(self.sess, full_path)
			return True, global_step
		else:
			return False, 0

	def test(self, noisy_total, eval_data, ckpt_dir, save_dir):
		"""Test Temp3-CNN"""
		# init variables
		tf.global_variables_initializer().run()
		assert len(eval_data) != 0, 'No testing data!'
		load_model_status, global_step = self.load(ckpt_dir)
		assert load_model_status == True, '[!] Load weights FAILED...'
		print(" [*] Load weights SUCCESS...")
		psnr_sum = 0
		ind = np.multiply(range(len(eval_data)/3),3)

		start = time.time()
		for idx in range(len(eval_data)):
			if idx>(len(eval_data)-3):
				test = cv2.imread(eval_data[idx])
				noisy = cv2.imread(noisy_total[idx])
				psnr = cal_psnr(test,noisy)
				psnr_sum += psnr
				print("img%d PSNR: %.2f %.2f" % (idx, psnr, psnr))
				_, path = eval_data[idx].split('al')
#				cv2.imwrite(('./data/denoised/%s' %path), noisy)
				print('last frames')
				continue
			test = cv2.imread(eval_data[idx])
			test1 = cv2.imread(eval_data[idx+1])
			test2 = cv2.imread(eval_data[idx+2])
			noisy = cv2.imread(noisy_total[idx])
			noisy1 = cv2.imread(noisy_total[idx+1])
			noisy2 = cv2.imread(noisy_total[idx+2])
			if idx==0:
				psnr = cal_psnr(test,noisy)
				psnr_sum += psnr
				_, path = eval_data[idx].split('al')
#				cv2.imwrite(('./data/denoised/%s' %path), noisy)

			test = test.astype(np.float32) / 255.0
			test1 = test1.astype(np.float32) / 255.0
			test2 = test2.astype(np.float32) / 255.0
			noisy = noisy.astype(np.float32) / 255.0
			noisy1 = noisy1.astype(np.float32) / 255.0
			noisy2 = noisy2.astype(np.float32) / 255.0

			orig = np.zeros((1,test.shape[0],test.shape[1],9))
			noisyin = np.zeros((1,test.shape[0],test.shape[1],9))

			orig[0,:,:,0:3] = test
			orig[0,:,:,3:6] = test1
			orig[0,:,:,6:] = test2
			noisyin[0,:,:,0:3] = noisy
			noisyin[0,:,:,3:6] = noisy1
			noisyin[0,:,:,6:] = noisy2

			output_clean_image= self.sess.run(
				[self.Y],feed_dict={self.Y_:orig,self.X:noisyin,self.is_training: False})

			out = np.asarray(output_clean_image)
#			cv2.imshow('',out[0,0])
#			cv2.waitKey(0)
			# calculate PSNR
			psnr = psnr_scaled(test1, out[0,0])
			psnr1 = psnr_scaled(test1, noisy1)
			print("img%d PSNR: %.2f %.2f" % (idx+1, psnr, psnr1))
			psnr_sum += psnr
			_, path = eval_data[idx+1].split('al')
#			cv2.imwrite(('./data/denoised/%s' %path), out[0,0]*255)
#			cv2.imwrite(('./data/original_frames/%s' %path), test1*255)
		avg_psnr = psnr_sum / len(eval_data)
		print("--- Average PSNR %.2f ---" % avg_psnr)
		print("--- Elapsed time: %.4f" %(time.time()-start))

def cal_psnr(im1, im2): # PSNR function for 0-255 values
	mse = ((im1.astype(np.float) - im2.astype(np.float)) ** 2).mean()
	psnr = 10 * np.log10(255 ** 2 / mse)
	return psnr

def psnr_scaled(im1, im2): # PSNR function for 0-1 values
	mse = ((im1 - im2) ** 2).mean()
	mse = mse * (255 ** 2)
	psnr = 10 * np.log10(255 **2 / mse)
	return psnr


text:
import numpy as np
import torch
import sys

# class GaussianKernelSpec(object):
#   def __init__(self, sigma):
#     self.sigma = sigma

class GaussianKernel(object):
  def __init__(self, sigma):
    self.sigma = sigma
    self.dist_func = torch.nn.PairwiseDistance(p=2)

  def get_kernel_matrix(self, X1, X2, quantizer1=None, quantizer2=None, dtype="float"):
    '''
    the input value has shape [n_sample, n_dim]
    quantizer is dummy here
    dtype only works for numpy input for X1 X2
    '''
    if isinstance(X1, np.ndarray) and isinstance(X2, np.ndarray):
      n_sample_X1 = X1.shape[0]
      norms_X1 = np.linalg.norm(X1, axis=1).reshape(n_sample_X1, 1)
      n_sample_X2 = X2.shape[0]
      norms_X2 = np.linalg.norm(X2, axis=1).reshape(n_sample_X2, 1)
      cross = np.dot(X1, X2.T)
      # print("using sigma ", self.sigma)
      kernel = np.exp(-0.5 / float(self.sigma)**2 \
        * (np.tile(norms_X1**2, (1, n_sample_X2) ) + np.tile( (norms_X2.T)**2, (n_sample_X1, 1) ) \
        -2 * cross) )
      if dtype == "float":
          return torch.Tensor(kernel).float()
      else:
          return torch.Tensor(kernel).double()
    else:
      ## to prevent memory explosion on GPU, we can do the following operations on CPU and move results
      ## back to GPU
      #is_cuda_tensor = X1.is_cuda
      #if is_cuda_tensor and use_cpu_comp:
      #    X1 = X1.cpu()
      #    X2 = X2.cpu()
      norms_X1 = (X1**2).sum(1).view(-1, 1)
      norms_X2 = (X2**2).sum(1).view(-1, 1)
      norms_X1 = norms_X1.repeat(1, int(X2.size(0) ) )
      norms_X2 = torch.transpose(norms_X2.repeat(1, int(X1.size(0) ) ), 0, 1)
      cross = torch.mm(X1, torch.transpose(X2, 0, 1) )
      kernel = torch.exp(-0.5 / float(self.sigma)**2 * (norms_X1 + norms_X2 - 2* cross) )
      #if is_cuda_tensor and use_cpu_comp:
      #    return kernel.cuda()
      #else:
      #    return kernel
      return kernel

  def torch(self, cuda=False):
    '''
    adapt the interface to the model launching wrapper
    '''
    pass

  def cpu(self):
    '''
    adapt the interface when switch parameter of some kernels back to cpu mode
    '''
    pass


text:
# -*- coding: utf-8 -*-

from psychopy import visual, event, core, monitors
from psychopy.hardware import joystick
import os, random
import numpy as np
from datetime import date

# Subject profile
# username = input("Please enter your name:").upper()
# print('Hi %s, welcome to our experiment!' % username)
# today = date.today()
# filename = ('%s_%s.txt' % (today, username))

# count = 0
# while os.path.isfile(filename):
#     count += 1
#     filename = ('%s_%s_%d.txt' % (today, username, count))

# f = open(filename, 'w')
# f.write('ohhhh heyyyy')
# f.close()

# print(filename)
# Make screen profile ----
widthPix = 2560 # screen width in px
heightPix = 2440 # screen height in px
monitorwidth = 60 # monitor width in cm
viewdist = 60 # viewing distance in cm
monitorname = 'ProArt27'
scrn = 0 # 0 to use main screen, 1 to use external screen
mon = monitors.Monitor(monitorname, width=monitorwidth, distance=viewdist)
mon.setSizePix((widthPix, heightPix))
mon.save()

# Make initials
experiment_timer = core.Clock()
experiment_timer.reset()
MAX_DURATION = 7
previous_resp = []
# Open window (for joystick)
my_win = visual.Window(size=(400, 400), pos=(0,0), monitor = mon, units = 'pix',
                       screen = 1)

# Get mouse and joystick
joystick.backend = 'pyglet'
nJoys = joystick.getNumJoysticks() # Check if I have any joysticks
id = 0 # Use the first one as input
joy = joystick.Joystick(id) # ID has to be nJoys - 1
# - Mouse setting
mouse = event.Mouse(visible = True, win = my_win)
mouse.clickReset() # Reset to its initials

# start
while experiment_timer.getTime() < MAX_DURATION:
    my_win.flip()
    a,b,c,d = getAnything(mouse, joy)

    if previous_resp != [a,b,c,d]:
        print(a,b,c,d)
        print(type(a), type(b), type(c), type(d))
    previous_resp = [a,b,c,d]

my_win.close()

max_stars_repo_name: EJ-Chang/OSD-Study

max_issues_repo_name: EJ-Chang/OSD-Study

max_forks_repo_name: EJ-Chang/OSD-Study


text:
import argparse
from cv2 import *
import numpy as np
arg=argparse.ArgumentParser()
arg.add_argument("--path",help="INPUT VIDEO PATH ",default=0)
args=vars(arg.parse_args())

cap=cv2.VideoCapture(args['path'])
while (True):
    check,frame=cap.read()
    gray=cvtColor(frame,COLOR_RGB2GRAY)
    img_brg=cvtColor(frame,COLOR_RGB2BGR)
    img_hsv=cvtColor(frame,COLOR_BGR2HSV)
    imshow("displaying GRAY IMAGE -->",gray)
    imshow("displaying IMAGE -->",frame)
    imshow("displaying BRG IMAGE -->", img_brg)
    imshow("displaying HSV IMAGE -->", img_hsv)
    if waitKey(100) & 0xFF==ord('q'):
        break;

cap.release()
destroyAllWindows()


text:
import h5py
import pickle
import time

import pandas as pd
import numpy as np

from tqdm import tqdm

from konlpy.tag import Okt

def explore(n):

    # Read columns of specified raw data that we would explore

    f = h5py.File('raw_data/train.chunk.0{}'.format(n), 'r')['train']
    keys = ['brand', 'maker', 'model', 'product', 'price']

    data = pd.DataFrame(columns=keys)

    for key in tqdm(keys):
        data[key] = f[key][:]

    #  Decode text columns to `utf-8`

    for key in tqdm(keys):
        if key != 'price':
            data[key] = data[key].str.decode('utf-8')

    #  Split texts from columns `brand` and `maker`

    data['brand'] = data['brand'].str.split('[^가-힣a-zA-Z]')
    data['maker'] = data['maker'].str.split('[^가-힣a-zA-Z]')

    # Morph texts from columns `model` and `product`

    data['model'] = data['model'].apply(Okt().morphs)
    data['product'] = data['product'].apply(Okt().morphs)

    # Replace `-1` with `np.nan` in `price` column

    data.loc[data[data['price']==-1].index, 'price'] = np.nan

    # Save to pickle format

    median_price = pd.Series(data['price'].median())

    maker_counts = pd.Series(np.concatenate(data['maker'])).value_counts()
    brand_counts = pd.Series(np.concatenate(data['brand'])).value_counts()
    model_counts = pd.Series(np.concatenate(data['model'])).value_counts()
    product_counts = pd.Series(np.concatenate(data['product'])).value_counts()

    total_counts = maker_counts.add(brand_counts, fill_value=0).add(model_counts, fill_value=0).add(product_counts, fill_value=0).sort_values(ascending=False)

    print(n)
    print(total_counts.head())
    print(median_price)

    with open('info/info-{}.pkl'.format(n), 'wb') as handle:
        pickle.dump((total_counts, median_price), handle)

for n in range(1,10):
    now = time.time()
    explore(n)
    print('TIME SPENT FOR ONE CHUNK: {}s'.format(int(time.time()-now)))


text:
#!/usr/bin/env python
# coding: utf-8
from __future__ import absolute_import
import numpy as np
import matplotlib.pyplot as plt
import matplotlib.backends.backend_pdf as backend_pdf

# Local imports
from ising.utils.file_utils import projfile, eigenfile, load_sparse_csr
import ising.utils.operator_utils as ou
import ising.utils.calculation_utils as cu
import ising.utils.ising_plot_utils as ipu
import ising.models.base as models
from examples.thermalization.params import Ls, MODELS, figLargeOpsPath

def plot_opstring(L, Op, k, evals, evecs):
    op_eev_mid, op_ev_mc, sigmaOp = cu.op_eev_fluct(L, Op, evals, evecs,
                                                    deltaE=.025*L,
                                                    spectrum_center=1/2,
                                                    spectrum_width=20)

    fig = plt.figure()
    plt.plot(op_eev_mid, '.', label='Eigenstate')
    plt.plot(op_ev_mc, label='Microcanonical')
    plt.title(r'Comparison to Microcanonical, $L=%d$' % L, fontsize=16)
    plt.title(r'$L=%d,~\Delta E=0.025L,$' % L
              + r'~$\mathcal{O}=(S^y)^{\otimes %d}\$' % k,
              fontsize=16)
    plt.ylim(-0.5, 0.5)
    plt.legend(fontsize=16)
    plt.xticks(fontsize=16)
    plt.yticks(fontsize=16)
    plt.tight_layout()

    return fig, op_eev_mid, op_ev_mc, sigmaOp

def plot_operators():
    pass

for model in MODELS:
    sigmaOp_vec = []
    for L in Ls:
        # ----------------------
        # Setting up model
        # ----------------------
        print("##################################\n"
              + "Calculations for model: "+str(model) + "; L="+str(L)
              + "\n##################################", flush=True)

        model_params, symmetries = models.default_params(model, L=L,
                                                         integrable=False)
        # Retrieving stored information
        print("    # ==============================\n"
              + "    # Retrieving stored information\n"
              + "    # ==============================", flush=True)
        proj_dict = load_sparse_csr(projfile(L, S=1/2, **symmetries))
        eigen_dict = np.load(eigenfile(model, **model_params),
                             allow_pickle=True)

        # Getting eigensystems
        print("        # --------------------------\n"
              + "        # Getting eigensystems\n"
              + "        # --------------------------\n", flush=True)
        sub_evals = eigen_dict['subspace evals'][()]
        sub_evecs = eigen_dict['subspace evecs'][()]
        all_evals, all_evecs = cu.esys_from_sub_dicts(proj_dict, eigen_dict)

        # Setting up operators
        sigma_list = ou.gen_s0sxsysz(L)[1:]
        # DEBUG: Don't need to do all of this for right now
        central_ops = ou.connected_central_ops(L, sigma_list, verbose=0)

        # Setting up lists to hold operators
        op_eevs_ys = []
        op_ev_mc_ys = []
        op_flucts_ys = []

        all_op_eevs = []
        all_op_ev_mcs = []
        all_op_flucts = []

        # -----------------
        # Finding and Plotting EEVs
        # -----------------
        # Setting up pdf with many plots
        ystring_file = backend_pdf.PdfPages(figLargeOpsPath+model+'_ystringev'
                                            + 'fixed_{}sites.pdf'.format(L))
        for k in range(L):
            y_string = central_ops[k]['y'*k]
            fig, eev, ev_mc, sigmaOp = ipu.plot_opstring(L, y_string, k,
                                                         all_evals, all_evecs)
            ystring_file.savefig(fig)

            op_eevs_ys.append(eev)
            op_ev_mc_ys.append(ev_mc)
            op_flucts_ys.append(sigmaOp)

            """
            all_op_eevs_k = []
            all_op_ev_mcs_k = []
            all_op_flucts_k = []

            # Enumerating over connected central operators
            for op_string in central_ops[k]:
                _, eevs = cu.op_ev(op_string, evecs)
                _, ev_mc, fluct = cu.op_eev_fluct(L, op_string, evals, evecs,
                                                  deltaE=.025*L)

                all_op_eevs_k.append(eevs)
                all_op_ev_mcs.append(ev_mc)
                all_op_flucts_k.append(fluct)

                # If the operator is a string of pauli xs,
                # consider it individually as well
                if central_ops[k].key(op_string) == 'x'*k:
                    op_eevs_xs.append(eevs)
                    op_ev_mc_xs.append(ev_mc)
                    op_flucts_xs.append(fluct)

            all_op_eevs.append(all_op_eevs_k)
            all_op_ev_mcs.append(all_op_ev_mcs_k)
            all_op_flucts.append(all_op_flucts_k)
            """
        ystring_file.close()

        fig = plt.figure()
        for k in range(L):
            fig.scatter([k]*len(op_eevs_ys), op_eevs_ys[k])
        fig.plot(range(L), op_ev_mc_ys)
        fig.fill_between(x=range(L), y1=op_ev_mc_ys+op_flucts_ys,
                         y2=op_ev_mc_ys-op_flucts_ys)

        fig.savefig(figLargeOpsPath+model+'_ystringev_{}sites.pdf'.format(L))

        # Store these lists in a dictionary of some sort

        """
        #---------------------------------
        # Plotting EEVs
        #---------------------------------
        # Set up figures:
        fig_eev_1 = aesthetic_fig()
        fig_flucts_1 = aesthetic_fig()
        fig_eev_all = aesthetic_fig()
        fig_flucts_all = aesthetic_fig()

        # Plotting
        # Could potentially set this up beforehand instead of storing data

        # Plotting for single operator
        fig_eev_1[0].plot(range(L), op_eevs_1)
        fig_flucts_1[0].plot(range(L), op_flucts_1)

        # Plotting for central operators
        for k in range(L):
            fig_eev_all[0].scatter(k, all_op_eevs[k])
            fig_flucts_all[0].scatter(k, all_op_flucts[k])
        """

max_stars_repo_name: samcaf/Ising

max_issues_repo_name: samcaf/Ising

max_forks_repo_name: samcaf/Ising


text:
# Copyright 2016 The TensorFlow Authors. All Rights Reserved.
#
# Licensed under the Apache License, Version 2.0 (the "License");
# you may not use this file except in compliance with the License.
# You may obtain a copy of the License at
#
#     http://www.apache.org/licenses/LICENSE-2.0
#
# Unless required by applicable law or agreed to in writing, software
# distributed under the License is distributed on an "AS IS" BASIS,
# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
# See the License for the specific language governing permissions and
# limitations under the License.
# ==============================================================================
"""Tests for layer wrappers."""

from __future__ import absolute_import
from __future__ import division
from __future__ import print_function

import numpy as np

from tensorflow.contrib.keras.python import keras
from tensorflow.python.platform import test

class TimeDistributedTest(test.TestCase):

  def test_timedistributed_dense(self):
    # first, test with Dense layer
    with self.test_session():
      model = keras.models.Sequential()
      model.add(
          keras.layers.TimeDistributed(
              keras.layers.Dense(2), input_shape=(3, 4)))
      model.add(keras.layers.Activation('relu'))
      model.compile(optimizer='rmsprop', loss='mse')
      model.fit(
          np.random.random((10, 3, 4)),
          np.random.random((10, 3, 2)),
          epochs=1,
          batch_size=10)

      # test config
      model.get_config()

  def test_timedistributed_conv2d(self):
    # test with Conv2D
    with self.test_session():
      model = keras.models.Sequential()
      model.add(
          keras.layers.TimeDistributed(
              keras.layers.Conv2D(5, (2, 2), padding='same'),
              input_shape=(2, 4, 4, 3)))
      model.add(keras.layers.Activation('relu'))
      model.compile(optimizer='rmsprop', loss='mse')
      model.train_on_batch(
          np.random.random((1, 2, 4, 4, 3)), np.random.random((1, 2, 4, 4, 5)))

      model = keras.models.model_from_json(model.to_json())
      model.summary()

  def test_timedistributed_stacked(self):
    # test stacked layers
    with self.test_session():
      model = keras.models.Sequential()
      model.add(
          keras.layers.TimeDistributed(
              keras.layers.Dense(2), input_shape=(3, 4)))
      model.add(keras.layers.TimeDistributed(keras.layers.Dense(3)))
      model.add(keras.layers.Activation('relu'))
      model.compile(optimizer='rmsprop', loss='mse')

      model.fit(
          np.random.random((10, 3, 4)),
          np.random.random((10, 3, 3)),
          epochs=1,
          batch_size=10)

  def test_timedistributed_sequential(self):
    # test wrapping Sequential model
    with self.test_session():
      model = keras.models.Sequential()
      model.add(keras.layers.Dense(3, input_dim=2))
      outer_model = keras.models.Sequential()
      outer_model.add(keras.layers.TimeDistributed(model, input_shape=(3, 2)))
      outer_model.compile(optimizer='rmsprop', loss='mse')
      outer_model.fit(
          np.random.random((10, 3, 2)),
          np.random.random((10, 3, 3)),
          epochs=1,
          batch_size=10)

      # test with functional API
      x = keras.layers.Input(shape=(3, 2))
      y = keras.layers.TimeDistributed(model)(x)
      outer_model = keras.models.Model(x, y)
      outer_model.compile(optimizer='rmsprop', loss='mse')
      outer_model.fit(
          np.random.random((10, 3, 2)),
          np.random.random((10, 3, 3)),
          epochs=1,
          batch_size=10)

  def test_regularizers(self):
    with self.test_session():
      model = keras.models.Sequential()
      model.add(
          keras.layers.TimeDistributed(
              keras.layers.Dense(2, kernel_regularizer='l1'),
              input_shape=(3, 4)))
      model.add(keras.layers.Activation('relu'))
      model.compile(optimizer='rmsprop', loss='mse')
      self.assertEqual(len(model.losses), 1)

  def test_TimeDistributed_learning_phase(self):
    # test layers that need learning_phase to be set
    np.random.seed(1234)
    x = keras.layers.Input(shape=(3, 2))
    y = keras.layers.TimeDistributed(
        keras.layers.Dropout(.999))(x, training=True)
    model = keras.models.Model(x, y)
    y = model.predict(np.random.random((10, 3, 2)))
    self.assertAllClose(np.mean(y), 0., atol=1e-1, rtol=1e-1)

class BidirectionalTest(test.TestCase):

  def test_bidirectional(self):
    rnn = keras.layers.SimpleRNN
    samples = 2
    dim = 2
    timesteps = 2
    output_dim = 2
    with self.test_session():
      for mode in ['sum', 'concat']:
        x = np.random.random((samples, timesteps, dim))
        target_dim = 2 * output_dim if mode == 'concat' else output_dim
        y = np.random.random((samples, target_dim))

        # test with Sequential model
        model = keras.models.Sequential()
        model.add(
            keras.layers.Bidirectional(
                rnn(output_dim), merge_mode=mode, input_shape=(timesteps, dim)))
        model.compile(loss='mse', optimizer='sgd')
        model.fit(x, y, epochs=1, batch_size=1)

        # test config
        model.get_config()
        model = keras.models.model_from_json(model.to_json())
        model.summary()

  def test_bidirectional_stacked(self):
    # test stacked bidirectional layers
    rnn = keras.layers.SimpleRNN
    samples = 2
    dim = 2
    timesteps = 2
    output_dim = 2
    mode = 'sum'

    with self.test_session():
      x = np.random.random((samples, timesteps, dim))
      target_dim = 2 * output_dim if mode == 'concat' else output_dim
      y = np.random.random((samples, target_dim))

      model = keras.models.Sequential()
      model.add(
          keras.layers.Bidirectional(
              rnn(output_dim, return_sequences=True),
              merge_mode=mode,
              input_shape=(timesteps, dim)))
      model.add(keras.layers.Bidirectional(rnn(output_dim), merge_mode=mode))
      model.compile(loss='mse', optimizer='sgd')
      model.fit(x, y, epochs=1, batch_size=1)

      # test with functional API
      inputs = keras.layers.Input((timesteps, dim))
      output = keras.layers.Bidirectional(
          rnn(output_dim), merge_mode=mode)(inputs)
      model = keras.models.Model(inputs, output)
      model.compile(loss='mse', optimizer='sgd')
      model.fit(x, y, epochs=1, batch_size=1)

  def test_bidirectional_statefulness(self):
    # Bidirectional and stateful
    rnn = keras.layers.SimpleRNN
    samples = 2
    dim = 2
    timesteps = 2
    output_dim = 2
    mode = 'sum'

    with self.test_session():
      x = np.random.random((samples, timesteps, dim))
      target_dim = 2 * output_dim if mode == 'concat' else output_dim
      y = np.random.random((samples, target_dim))

      inputs = keras.layers.Input(batch_shape=(1, timesteps, dim))
      output = keras.layers.Bidirectional(
          rnn(output_dim, stateful=True), merge_mode=mode)(inputs)
      model = keras.models.Model(inputs, output)
      model.compile(loss='mse', optimizer='sgd')
      model.fit(x, y, epochs=1, batch_size=1)

if __name__ == '__main__':
  test.main()

max_stars_repo_name: DEVESHTARASIA/tensorflow

max_issues_repo_name: DEVESHTARASIA/tensorflow

max_forks_repo_name: DEVESHTARASIA/tensorflow


text:
import matplotlib.pyplot as plt
import numpy as np
from scipy.io import savemat

# Hardcoded binary image for illustration of a heart
sample_x = np.array([[1, 0, 1, 0, 1],
		[0, 1, 0, 1, 0],
		[0, 1, 1, 1, 0],
		[1, 0, 1, 0, 1],
		[1, 1, 0, 1, 1]])

x, y = sample_x.shape
string_list = []
for i in range(x):
	for j in range(y):
		string_list.append("{},{}".format(i, j))
indices_A = np.array(string_list, dtype = 'object').reshape(x, y)

# Random alpha scaling to prevent trivial 1, 0 problem
# import random
# alpha_scale = random.sample(range(1, 8), y)
# for i, alpha in enumerate(alpha_scale):
# 	sample_x[:, i] = alpha * sample_x[:, i]

# Alternatively hardcoded intensity scaling for reproducibility
sample_x[:, 1] =  5*sample_x[:, 1]
sample_x[:, 2] =  2*sample_x[:, 2]
sample_x[:, 3] =  3*sample_x[:, 3]
sample_x[:, 4] =  6*sample_x[:, 4]
sample_x[:, 0] =  4*sample_x[:, 0]

# Obtain actual vector x used in model
x_flatten = sample_x.flatten()

# Visualize the image
print(sample_x)
plt.figure()
plt.imshow(sample_x)
plt.title('Sample Image')
plt.show()
plt.close()

# Computing rays along all possible diagonals and their respective indices ('row', 'column')
def sample_diag(sample_x):
	diags = [sample_x[::-1, :].diagonal(i) for i in range(-sample_x.shape[0]+1, sample_x.shape[1])]
	diags.extend(sample_x.diagonal(i) for i in range(sample_x.shape[1]-1, -sample_x.shape[0], -1))
	return [n.tolist() for n in diags]
all_diags = sample_diag(sample_x)
diag_indic = sample_diag(indices_A)

# Computing rays along all columns
column_arr = []
column_indic_global = []
for j in range(y):
	columns = [row[j] for row in sample_x]
	column_arr.append(columns)
	column_indic = []
	for r_cnt in range(x):
		column_indic.append("{},{}".format(r_cnt, j))
	column_indic_global.append(column_indic)

# Computing rays along all rows
rows_arr = []
row_indic_global = []
for j in range(x):
	rows = [row for row in sample_x[j, :]]
	rows_arr.append(rows)
	row_indic = []
	for c_cnt in range(y):
		row_indic.append("{},{}".format(j, c_cnt))
	row_indic_global.append(row_indic)

# Joining diagonals, columns and row indices and values into a single list
all_indices = [*diag_indic, *column_indic_global, *row_indic_global]
all_values_combo = [*all_diags, *column_arr, *rows_arr]

# Creating empty mixing matrix A for model
row_dim, col_dim = np.shape(sample_x)
A_mixing_matrix = np.zeros((len(all_indices), row_dim, col_dim))
for obs_n, ind in enumerate(all_indices):
	for sub_i in ind:
		local_x, local_y = sub_i.split(',')
		A_mixing_matrix[obs_n, int(local_x), int(local_y)] = 1
A_final = A_mixing_matrix.reshape((obs_n+1, row_dim*col_dim))

# Pure observations are used for lambda to sample from poisson distribution
pure_obser = np.sum(A_final*x_flatten, axis=1)
poisson_obser = np.random.poisson(pure_obser, len(all_indices))

# Saving to mat files for use in MATLAB
export_dict = {"A": A_final, "y":poisson_obser, "x": x_flatten}

savemat("simulated_heart.mat", export_dict)
print("Image vector x is of shape {}".format(len(x_flatten)))
print("Mixing matrix A is of shape {}".format(np.shape(A_final)))
print("Obervation matrix Y is of shape {}".format(np.shape(poisson_obser)))
print("X:", x_flatten)
print("A:", A_final)
print("Y:", poisson_obser)
print("Y without poisson", pure_obser)

reason: import numpy,from scipy


text:
""" Testing reading DICOM files

"""

import numpy as np

from .. import dicomreaders as didr

from nibabel.pydicom_compat import dicom_test, pydicom

from .test_dicomwrappers import (EXPECTED_AFFINE,
                                 EXPECTED_PARAMS,
                                 IO_DATA_PATH,
                                 DATA)

from nose.tools import (assert_true, assert_false, assert_equal, assert_raises)

from numpy.testing import assert_array_equal, assert_array_almost_equal

@dicom_test
def test_read_dwi():
    img = didr.mosaic_to_nii(DATA)
    arr = img.get_data()
    assert_equal(arr.shape, (128, 128, 48))
    assert_array_almost_equal(img.affine, EXPECTED_AFFINE)

@dicom_test
def test_read_dwis():
    data, aff, bs, gs = didr.read_mosaic_dwi_dir(IO_DATA_PATH,
                                                 'siemens_dwi_*.dcm.gz')
    assert_equal(data.ndim, 4)
    assert_array_almost_equal(aff, EXPECTED_AFFINE)
    assert_array_almost_equal(bs, (0, EXPECTED_PARAMS[0]))
    assert_array_almost_equal(gs, (np.zeros((3,)), EXPECTED_PARAMS[1]))
    assert_raises(IOError, didr.read_mosaic_dwi_dir, 'improbable')

@dicom_test
def test_passing_kwds():
    # Check that we correctly pass keywords to dicom
    dwi_glob = 'siemens_dwi_*.dcm.gz'
    csa_glob = 'csa*.bin'
    for func in (didr.read_mosaic_dwi_dir, didr.read_mosaic_dir):
        data, aff, bs, gs = func(IO_DATA_PATH, dwi_glob)
        # This should not raise an error
        data2, aff2, bs2, gs2 = func(
            IO_DATA_PATH,
            dwi_glob,
            dicom_kwargs=dict(force=True))
        assert_array_equal(data, data2)
        # This should raise an error in pydicom.dicomio.read_file
        assert_raises(TypeError,
                      func,
                      IO_DATA_PATH,
                      dwi_glob,
                      dicom_kwargs=dict(not_a_parameter=True))
        # These are invalid dicoms, so will raise an error unless force=True
        assert_raises(pydicom.filereader.InvalidDicomError,
                      func,
                      IO_DATA_PATH,
                      csa_glob)
        # But here, we catch the error because the dicoms are in the wrong
        # format
        assert_raises(didr.DicomReadError,
                      func,
                      IO_DATA_PATH,
                      csa_glob,
                      dicom_kwargs=dict(force=True))

max_stars_repo_name: Raniac/NEURO-LEARN

max_forks_repo_name: Raniac/NEURO-LEARN

reason: import numpy,from numpy


text:
import re
from collections import Counter
from dataclasses import dataclass
from enum import IntFlag
from functools import reduce
from operator import ior
from typing import Dict, List, Optional, Set, Tuple, Union
from functools import lru_cache
import networkx as nx

import utils

class NodeAttribute(IntFlag):
    REGULAR = 0
    TREATMENT = 1
    OUTCOME = 2
    ADJUSTED = 4
    UNOBSERVED = 8

    @classmethod
    def parse(cls, attr: str) -> "NodeAttribute":
        attr = attr.lower()
        if attr == "regular":
            return cls.REGULAR
        elif attr in ["treatment", "t"]:
            return cls.TREATMENT
        elif attr in ["outcome", "o"]:
            return cls.OUTCOME
        elif attr in ["adjusted", "a"]:
            return cls.ADJUSTED
        elif attr in ["unobserved", "u"]:
            return cls.UNOBSERVED
        else:
            raise ValueError(f"Can't parse '{attr}'")

    def __str__(self) -> str:
        if self == NodeAttribute.REGULAR:
            return ""
        else:
            return self.name.lower()

@dataclass
class CausalGraph:
    nodes: Set[str] = None
    edges: Set = None
    treatment: Optional[str] = None
    outcome: Optional[str] = None
    adjusted: Set[str] = None
    unobserved: Set[str] = None
    # compounds: Dict[Set[str]] = None
    _graph: nx.DiGraph = None

    def __post_init__(self):
        self._build_graph()

    def _build_graph(self):
        self._graph = nx.DiGraph()

        if self.nodes is not None:
            self._graph.add_nodes_from(self.nodes, observed=True, parent="all")
        if self.edges is not None:
            self._graph.add_edges_from(self.edges)

        # set node attributes
        node_attrs = {}
        if self.treatment is not None:
            node_attrs[self.treatment] = {"treatment": True}
        if self.outcome is not None:
            node_attrs[self.outcome] = {"outcome": True}
        if self.adjusted is not None:
            node_attrs.update({n: {"adjusted": True} for n in self.adjusted})
        if self.unobserved is not None:
            node_attrs.update({n: {"observed": False} for n in self.unobserved})
        nx.set_node_attributes(self._graph, node_attrs)
        self.update_paths()

    @property
    def graph(self):
        return self._graph

    def update_paths(self):
        causal_edges = [e for path in self.get_causal_paths(as_edge_list=True) for e in path]
        biasing_edges = [e for path in self.get_biasing_paths(as_edge_list=True) for e in path]

        edge_attrs = {e: {
            "causal": e in causal_edges,
            "biasing": e in biasing_edges
        } for e in self._graph.edges()}

        nx.set_edge_attributes(self._graph, edge_attrs)

    def get_node_attributes(self, node: str) -> NodeAttribute:
        if self.treatment == node:
            return NodeAttribute.TREATMENT
        elif self.outcome == node:
            return NodeAttribute.OUTCOME
        elif node in self.adjusted:
            return NodeAttribute.ADJUSTED
        elif node in self.unobserved:
            return NodeAttribute.UNOBSERVED
        else:
            return NodeAttribute.REGULAR

    def add_node(self, node_id: str) -> None:
        self.nodes.add(node_id)
        self._graph.add_node(node_id, observed=True)

    def add_edge(self, source: str, target: str) -> None:
        self.edges.add((source, target))
        self._graph.add_edge(source, target)

    def update_node(self, node: str, attr: NodeAttribute) -> str:
        new_node_attrs = None
        if attr == NodeAttribute.TREATMENT:
            new_node_attrs = {
                self.treatment: {"treatment": False},
                node: {"treatment": True}
            }
            self.treatment = node
        elif attr == NodeAttribute.OUTCOME:
            new_node_attrs = {
                self.outcome: {"outcome": False},
                node: {"outcome": True}
            }
            self.outcome = node
        elif attr == NodeAttribute.ADJUSTED:
            new_node_attrs = {node: {"adjusted": True}}
            self.adjusted.add(node)
        elif attr == NodeAttribute.UNOBSERVED:
            self.unobserved.add(node)
            new_node_attrs = {node: {"observed": False}}
        elif attr == NodeAttribute.REGULAR:
            # reset all other node attributes
            new_node_attrs = {node: {
                "treatment": False,
                "outcome": False,
                "adjusted": False,
                "observed": True,
            }}
            if node in self.adjusted:
                self.adjusted.remove(node)
            if node in self.unobserved:
                self.unobserved.remove(node)

        if new_node_attrs is not None:
            nx.set_node_attributes(self._graph, new_node_attrs)

    def set_node_position(self, node: str, x, y) -> None:
        node_attr = {"position": {"x": x, "y": y}}
        nx.set_node_attributes(self._graph, {node: node_attr})

    def delete_node(self, node_id: str) -> None:
        if self.treatment == node_id:
            self.treatment = None
        if self.outcome == node_id:
            self.outcome = None
        if node_id in self.adjusted:
            self.adjusted.remove(node_id)
        if node_id in self.unobserved:
            self.unobserved.remove(node_id)
        self.nodes.remove(node_id)
        self.edges = {(s, t) for s, t in self.edges if s != node_id and t != node_id}
        self._graph.remove_node(node_id)

    def delete_edge(self, source: str, target: str) -> None:
        self.edges.remove((source, target))
        self._graph.remove_edge(source, target)

    def as_string(self) -> str:
        edges = [f"{s}[{str(self.get_node_attributes(s))}]->{t}[{str(self.get_node_attributes(t))}]" for s, t in
                 list(self.edges)]
        return "\n".join(edges)

    def get_causal_paths(self, as_edge_list: bool = False) -> List[List[Tuple]]:
        if self.treatment is None or self.outcome is None:
            return []

        if as_edge_list:
            paths = nx.all_simple_edge_paths(self._graph, self.treatment, self.outcome)
        else:
            paths = nx.all_simple_paths(self._graph, self.treatment, self.outcome)
        return list(paths)

    def get_biasing_paths(self, as_edge_list: bool = True) -> List[List[Tuple]]:
        if self.treatment is None or self.outcome is None:
            return []
        backdoor_paths = self.get_backdoor_paths(self.treatment, self.outcome, as_edge_list=True)
        biasing_paths = [path for path in backdoor_paths if not self.is_path_blocked(path, self.adjusted)]
        if not as_edge_list:
            biasing_paths = [utils.edge_path_to_node_path(p) for p in biasing_paths]
        return biasing_paths

    def get_backdoor_paths(self, source: str, target: str, as_edge_list: bool = False):
        undirected_graph = self._graph.to_undirected()
        # by definition a backdoor path is any path to the target node, which starts
        # with an edge pointing to the source node (i.e. back door)
        backdoor_paths = [
            pth
            for pth in nx.all_simple_paths(undirected_graph, source=source, target=target)
            if self._graph.has_edge(pth[1], pth[0])]

        if as_edge_list:
            backdoor_edge_paths = []
            # the identified paths are undirected; the keep consistency, the edges
            # on the path have to be in the right direction
            for path in backdoor_paths:
                fixed_path = []
                for s, t in zip(path, path[1:]):
                    fixed_path.append((s, t) if self._graph.has_edge(s, t) else (t, s))
                backdoor_edge_paths.append(fixed_path)
            backdoor_paths = backdoor_edge_paths
        return backdoor_paths

    def is_path_blocked(self, path: List[Tuple], conditioning_set: Optional[Set] = None) -> bool:
        # rule 1 (unconditional separation): is there no collider
        target_counts = Counter([target for src, target in path])
        path_nodes = set([node for edge in path for node in edge])
        colliders = {node for node, c in target_counts.items() if c > 1}
        if conditioning_set is None or len(conditioning_set) == 0:
            return len(colliders) > 0

        # rule 2 (blocking by conditioning): is a variable on the path conditioned on
        # only variables on the path are valid blockers
        blocking_nodes = conditioning_set.union(colliders).intersection(path_nodes)

        # rule 3 (conditioning on colliders): paths blocked by colliders
        # which are conditioned on (or its descendants) are 'open'
        collider_descendants = {n: nx.descendants(self._graph, n) for n in colliders}
        descendants = reduce(set.union, collider_descendants.values(), set())

        defused_colliders = colliders.union(descendants).intersection(conditioning_set)
        blocking_nodes = blocking_nodes.difference(defused_colliders)

        return len(blocking_nodes) > 0

    def get_post_treatment_nodes(self):
        return nx.descendants(self._graph, self.treatment)

    def get_unobserved_nodes(self):
        return self.unobserved

def normalize_node(node_name: str) -> Tuple[str, NodeAttribute]:
    attrs = NodeAttribute.REGULAR
    match = re.search(r"(.*)\[(.*)\]", node_name)
    if match:
        node_name, attributes = match.groups()
        try:
            attrs = reduce(ior, [NodeAttribute.parse(a) for a in attributes.split(",") if a], NodeAttribute.REGULAR)
        except ValueError as exc:
            print(f"Can't parse node attributes: '{attributes}': {exc}")
    return node_name, attrs

def parse_model_string(model_string: Union[str, List[str]]) -> CausalGraph:
    treatment = None
    outcome = None
    adjusted = set()
    unobserved = set()
    nodes = set()
    edges = set()
    compounds = {}

    lines = model_string.split("\n") if isinstance(model_string, str) else model_string
    for line in lines:
        is_compound = False
        if "->" in line:
            source, target, *_ = [n.strip() for n in line.split("->")]
        elif "<-" in line:
            target, source, *_ = [n.strip() for n in line.split("<-")]
        elif "∈" in line:
            source, target, *_ = [n.strip() for n in line.split("∈")]
            is_compound = True
        else:
            raise ValueError(f"No valid edge direction specified in '{line}'. Expecting either '->' or '<-'")

        if len(_) > 0:
            raise NotImplementedError("Parsing of more than one edge per line are not yet supported!")

        source = normalize_node(source)
        target = normalize_node(target)

        if is_compound:
            parent = compounds.get(target, set())
            parent.add(source)
            compounds[target] = parent
        else:
            edges.add((source[0], target[0]))

        for node, node_attr in [source, target]:
            nodes.add(node)
            if node_attr == NodeAttribute.TREATMENT:
                if treatment is not None and treatment != node:
                    raise ValueError(
                        f"There is already a defined treatment: '{treatment}''. It will be overwritten with {node}")
                treatment = node
            if node_attr == NodeAttribute.OUTCOME:
                if outcome is not None and outcome != node:
                    raise ValueError(
                        f"There is already a defined outcome: '{outcome}'. It will be overwritten with {node}")
                outcome = node
            if node_attr == NodeAttribute.UNOBSERVED:
                unobserved.add(node)
            if node_attr == NodeAttribute.ADJUSTED:
                adjusted.add(node)

    model = CausalGraph(nodes, edges, treatment, outcome, adjusted, unobserved)
    return model

max_stars_repo_name: Magier/Aetia

max_issues_repo_name: Magier/Aetia

max_forks_repo_name: Magier/Aetia

reason: import networkx


text:
import dash
from dash.dependencies import Input, Output, State

import dash_core_components as dcc
import dash_html_components as html
import plotly.express as px
import plotly.graph_objects as go
import dash_bootstrap_components as dbc
import numpy as np
from utils.helper import *
import logging
from random import randrange
import pandas as pd
import dash_table
from datetime import datetime as dt
from datetime import date, timedelta

external_stylesheets = ['https://codepen.io/chriddyp/pen/bWLwgP.css']

# app = dash.Dash(__name__, external_stylesheets=external_stylesheets)
app = dash.Dash(__name__,
                external_stylesheets=[dbc.themes.BOOTSTRAP]
                )

app.layout = html.Div(children=[
    html.Div(
        id="sidebar",
        className="sidebar",
        children=[
            html.H5(
                style={"color": "white"},
                children=["SALESBOT"]
            ),
            html.A(id="a1_title_sidebar", className="title_sidebar", children=["Overall Performance"]),
            html.A(id="a2_title_sidebar", className="title_sidebar", children=["Bot performance by Outcomes"]),
            html.A(id="a3_title_sidebar", className="title_sidebar", children=["Bot performance by Use cases"]),
            html.A(id="a4_title_sidebar", className="title_sidebar", children=["Outcome of each Use case"]),
        ]
    ),

    html.Div(
        id="page_content_and_graph",
        children=[
            html.Div(
                id="part_1_title",
                children=["Overall Performance"],
            ),
            html.Div(
                id="part_1_subtitle_1",
                className="part_1_subtitle",
                children=["Statistics by last week: abc"],
            ),
            html.Div(
                className="basic-metrics",
                children=[
                    html.Div(
                        className="col-md-4",
                        children=[
                            html.Div(
                                id="no_conversations_lastweek",
                                className="sub_basic_metrics",
                                children=[
                                    html.Div(
                                        className="col-md-12",
                                        children=[
                                            html.Img(src="assets/icon/conversation_icon.png",
                                                     className="sub_basic_metrics_img"),
                                            html.P("Conversations", style={"display": "inline", "paddingLeft": "19px",
                                                                           "fontSize": "19px"}),
                                        ]
                                    )
                                ]
                            ),
                        ]
                    ),
                    html.Div(
                        className="col-md-4",
                        children=[
                            html.Div(
                                id="no_users_lastweek",
                                className="sub_basic_metrics",
                                children=[
                                    html.Div(
                                        className="col-md-12",
                                        children=[
                                            html.Img(src="assets/icon/user_icon.png",
                                                     className="sub_basic_metrics_img"),
                                            html.P("Users", style={"display": "inline", "paddingLeft": "19px",
                                                                   "fontSize": "19px"}),
                                        ]
                                    )
                                ]
                            ),
                        ]
                    ),
                    html.Div(
                        className="col-md-4",
                        children=[
                            html.Div(
                                id="success_rate_lastweek",
                                className="sub_basic_metrics",
                                children=[
                                    html.Div(
                                        className="col-md-12",
                                        children=[
                                            html.Img(src="assets/icon/success_rate_icon.png",
                                                     className="sub_basic_metrics_img"),
                                            html.P("Success Rate", style={"display": "inline", "paddingLeft": "19px",
                                                                          "fontSize": "19px"}),
                                        ]
                                    )
                                ]
                            ),
                        ]
                    ),
                ]
            ),
            html.Div(
                id="part_1_subtitle_2",
                className="part_1_subtitle",
                children=["Statistics by month"],
            ),
            html.Div(
                className="basic-metrics",
                children=[
                    html.Div(
                        className="col-md-4",
                        children=[
                            html.Div(
                                id="no_conversations_by_month",
                                className="sub_basic_metrics",
                                children=[
                                    html.Div(
                                        className="col-md-12",
                                        children=[
                                            # html.Img(src="assets/icon/conversation_icon.png", className="sub_basic_metrics_img"),
                                            html.P("Conversations by month",
                                                   style={"display": "inline", "paddingLeft": "19px",
                                                          "fontSize": "19px"}),
                                        ]
                                    )
                                ]
                            ),
                        ]
                    ),
                    html.Div(
                        className="col-md-4",
                        children=[
                            html.Div(
                                id="no_users_by_month",
                                className="sub_basic_metrics",
                                children=[
                                    html.Div(
                                        className="col-md-12",
                                        children=[
                                            html.P("Users by month", style={"display": "inline", "paddingLeft": "19px",
                                                                            "fontSize": "19px"}),
                                        ]
                                    )
                                ]
                            ),
                        ]
                    ),
                    html.Div(
                        className="col-md-4",
                        children=[
                            html.Div(
                                id="success_rate_by_month",
                                className="sub_basic_metrics",
                                children=[
                                    html.Div(
                                        className="col-md-12",
                                        children=[
                                            html.P("Success rate by month",
                                                   style={"display": "inline", "paddingLeft": "19px",
                                                          "fontSize": "19px"}),
                                        ]
                                    )
                                ]
                            ),
                        ]
                    ),
                ]
            ),
            html.Div(
                id="part_1_subtitle_3",
                className="part_1_subtitle",
                children=["Overall Performance by selected period"],
            ),
            html.Div(
                children=[
                    dcc.DatePickerSingle(
                        id='my_date_picker_start',
                        min_date_allowed=dt(2020, 1, 1),
                        max_date_allowed=dt(2021, 12, 31),
                        style={"borderRadius": "25px"},
                        placeholder='Start date',
                    ),
                    dcc.DatePickerSingle(
                        style={"paddingLeft": "20px", "borderRadius": "25px", "paddingRight": "20px"},
                        id='my_date_picker_end',
                        min_date_allowed=dt(2020, 1, 1),
                        max_date_allowed=dt(2021, 12, 31),
                        placeholder='End date',
                    ),
                    html.Button('Run', id='run-analytics', style={"borderRadius": "4px", "backgroundColor": "#448efc", "color": "white"}, ),
                ]
            )

        ],
    ),

])

if __name__ == '__main__':
    app.run_server(debug=True)