text:
import numpy as np
from astropy import units as u
from astropy import constants as const
from astrodendro import Dendrogram
from astropy.io import fits
from astropy.io.fits import getheader, getdata
import os

def define_asgn(image,dendrogram,label_out='',write=True,check_structures=False):
    '''
    PURPOSE: Writes and returns a 3D fits file with the index from the dendrogram as the BVALUE
        Required parameters:
            image: location of fits cube
            dendrogram: string/location of dendrogram hdf5 file
        Optional keywords:
            label_out: string to save assignment cube. automatically adds .asgn.fits.gz to end of string,
                       defaults to {current working directory}+'assignment_cube.asgn.fits.gz'
            write: set to False to prevent writing fits file
            check_structures: set to True to check if all of the structures in the dendrogram have
                              been written to the assignment cube (will take a long time for large dendrograms)
    '''
    if label_out == '':
        label_out = os.getcwd()+'/assignment_cube'
    d = Dendrogram.load_from(dendrogram)
    cube, hd3 = getdata(image, header=True)
    asgn = np.ones(cube.shape).astype(np.int)
    asgn *= -1
    def recursive_structures(trunks):
        for trunk in trunks:
            j = trunk.idx
            asgn[d[j].get_mask(shape=asgn.shape,subtree=False)] = j
            if len(trunk.children) > 0:
                recursive_structures(trunk.children)
    trunks = d.trunk
    recursive_structures(trunks)
    if check_structures:
        all_structures_present(asgn, d)
    if write:
        new_header = hd3
        new_header.set('BUNIT','Index')
        hdu = fits.PrimaryHDU(asgn)
        hdu.header = new_header
        hdu.header['DATAMIN'] = -1
        hdu.header['DATAMAX'] = len(d) - 1
        hdu.writeto(label_out+'.asgn.fits.gz')
        print('Wrote assignment cube to ', label_out)
    return asgn

def all_structures_present(asgn, dendrogram):
    '''
    PURPOSE: Check if all structures in the dendrogram are present in the assignment array
        Required parameters:
            asgn: assignment array
            dendrogram: Dendrogram object
    '''
    num = 0
    for i in range(len(dendrogram)):
        if len(np.where(asgn == i)[0]) > 0:
            num += 1
    if num != len(dendrogram):
        print('not all structures present in asgn cube')
        return False
    else:
        print('all structures present in asgn cube')
        return True

reason: import numpy,from astropy


text:
import collections
from typing import Tuple, List, Any, Dict, Optional

import gym
import numpy as np
import pandas as pd
from gym import spaces
from gym import utils

ITEM_METADATA_KEY = "item_metadata"

class RecSysEnv(gym.Env, utils.EzPickle):
    metadata = {"render.modes": []}
    reward_range = [0.0, 1.0]

    def __init__(
        self,
        dataset: pd.DataFrame,
        item_column: str,
        number_of_items: int,
        available_items_column: Optional[str] = None,
        item_metadata: Optional[Dict[str, np.ndarray]] = None,
    ):
        super().__init__()
        self._dataset = dataset.copy()
        self._item_metadata = item_metadata
        self._item_column = item_column
        self._available_items_column = available_items_column

        self._number_of_items = (
            number_of_items  # number_of_itemsdataset[item_column].max() + 1
        )

        if available_items_column:
            assert isinstance(
                self._dataset[available_items_column].values[0], collections.Sequence
            )
            available_items = np.zeros(
                (len(self._dataset), self._number_of_items), dtype=np.int8
            )
            non_zero_indices = [
                (i, available_item)
                for i, available_items in enumerate(
                    self._dataset[available_items_column]
                )
                for available_item in available_items
            ]
            i, j = zip(*non_zero_indices)
            available_items[i, j] = 1
            self._dataset[available_items_column] = list(iter(available_items))

        self._obs_dataset: List[dict] = self._dataset.drop(
            columns=[item_column]
        ).to_dict("records")

        self.action_space = spaces.Discrete(self._number_of_items)

        observation_space = {
            key: self._convert_value_to_space(key, value)
            for key, value in self._obs_dataset[0].items()
        }
        if item_metadata is not None:
            observation_space[ITEM_METADATA_KEY] = spaces.Dict(
                {
                    key: spaces.Box(
                        value.min(), value.max(), shape=value.shape, dtype=value.dtype
                    )
                    for key, value in item_metadata.items()
                }
            )
        self.observation_space = spaces.Dict(observation_space)
        self._current_index = 0

    def _convert_value_to_space(self, key: str, value: Any) -> spaces.Space:
        if isinstance(value, list):
            value = np.array(value)

        if key == self._available_items_column:
            return spaces.MultiBinary(self._number_of_items)
        if isinstance(value, int):
            return spaces.Discrete(self._dataset[key].max() + 1)
        elif isinstance(value, float):
            return spaces.Box(
                self._dataset[key].min(), self._dataset[key].max(), shape=(1,)
            )
        elif isinstance(value, np.ndarray):
            if issubclass(value.dtype.type, np.integer):
                return spaces.MultiDiscrete(
                    [np.max(self._dataset[key].max()) + 1] * len(value)
                )
            elif issubclass(value.dtype.type, np.floating):
                return spaces.Box(
                    self._dataset[key].min(),
                    self._dataset[key].max(),
                    shape=value.shape,
                )
        raise ValueError(
            "Unkown type in the observation space for {}:{}".format(key, value)
        )

    def _compute_stats(self, action: float) -> dict:
        # TODO: Choose which batch metrics to return
        return {}

    def _compute_reward(self, action: int) -> float:

        return float(
            self._dataset.iloc[self._current_index][self._item_column] == action
        )

    def _get_next_ob(self) -> dict:
        ob = self._obs_dataset[self._current_index].copy()
        if self._item_metadata is not None:
            ob[ITEM_METADATA_KEY] = self._item_metadata
        else:
            ob[ITEM_METADATA_KEY] = None
        return ob

    def step(self, action: int) -> Tuple[dict, float, bool, dict]:
        reward = self._compute_reward(action)
        info = self._compute_stats(action)

        self._current_index += 1

        done = (self._current_index + 1) == len(self._dataset)
        next_ob = self._get_next_ob() if not done else None

        return next_ob, reward, done, info

    def reset(self) -> dict:
        self._current_index = 0
        return self._get_next_ob()

    def render(self, mode="human"):
        pass

    def close(self):
        pass

ext: py

max_stars_repo_name: engsoares/mars-gym

max_issues_repo_name: kiminh/mars-gym

max_forks_repo_name: kiminh/mars-gym

reason: import numpy


text:
# -*- coding: utf-8 -*-
# Copyright(c) 2009 - present CNRS
# All rights reserved.

import argparse
import bloscpack as bp
from sklearn.svm import LinearSVC
from sklearn.metrics import accuracy_score,classification_report,confusion_matrix
from sklearn.model_selection import GridSearchCV as gsc
from sklearn.model_selection import RandomizedSearchCV as rsc
import scipy.sparse as sp
import numpy as np
from scipy.stats import expon
import matplotlib.pyplot as plt
from sklearn.metrics import plot_confusion_matrix

parser = argparse.ArgumentParser(
    description='SVM classifier using pre-trained BERT sentence representations  ')

parser.add_argument('--input', type=str, required=True, help='Path to the dataset')
parser.add_argument('--output', type=str, required=True,help="directory to save the output files")
parser.add_argument('--lang', type=str, required=True,help="fr or zh")
args = parser.parse_args()

X_train = bp.unpack_ndarray_from_file(args.input+'/train_X.blp')
Y_train = bp.unpack_ndarray_from_file(args.input+'/train_Y.blp')

X_test = bp.unpack_ndarray_from_file(args.input+'/test_X.blp')
Y_test = bp.unpack_ndarray_from_file(args.input+'/test_Y.blp')

# train the model on training set
model = LinearSVC(class_weight="balanced",random_state=42)
# defining parameter range
params = {"C": [0.0001,0.001,0.01,0.02,0.03,0.04,0.05,0.1]}
tense_clf = gsc(model,params,n_jobs=5,refit = True)

#fitting the model for grid search
tense_clf.fit(X_train,Y_train)

# print best parameter after tuning
print(tense_clf.best_params_)
# print how our model looks after hyper-parameter tuning
print(tense_clf.best_estimator_)

y_pred= tense_clf.predict(X_test)
print(classification_report(Y_test,y_pred))
#plot_confusion_matrix(model,X_test,Y_test,display_labels=['Past','Fut','Pres'],cmap=plt.cm.Blues,normalize='true')

reason: import numpy,import scipy,from scipy


text:
"""Flicker Plotting Functions

These functions generate commonly-used flicker graphics.

The functions are:

    * ieee_par_1789_graph - Plots the IEEE PAR 1789 logarithmic graph
    * waveform_graph - Plots the time-domain flicker waveform
    * standards_color - Returns colors for decorating standards result labels
"""

import math
import numpy as np
import matplotlib.pyplot as plt
import matplotlib.colors
import itertools
from matplotlib.ticker import PercentFormatter, ScalarFormatter
from .utils import bool_to_pass_fail
from sklearn.preprocessing import minmax_scale
from pylab import text

def ieee_par_1789_graph(
        data, figsize:tuple=(8,4), filename:str=None, showred:bool=True, showyellow:bool=True,
        noriskcolor:bool=True, max_freq:int=3000, min_pct:float=0.1, suppress:bool=False
    ):
    """Plots the IEEE PAR 1789 logarithmic graph

    Parameters
    ----------
    data : list
        The data points as a list of tuples as [(freq, mod, 'Name')], where mod is <= 1
    figsize : tuple
        The (width,height) of the plotted figure
    filename : str or None
        If specified, will save plot as the specified filename, e.g.: filename='../out/this_graph.png'
    showred : bool
        Whether to show the unsafe region in red
    showyellow : bool
        Whether to show the low-risk region in yellow
    noriskcolor : bool
        If False, no risk region will show in gray. If True, will show in green
    max_freq : int
        The maximum frequency, in kHz
    min_pct : float
        The minimum percent to display
    suppress : bool
        If True, the plot will not be shown
    """

    # count minimum percent decimals and recompute
    if min_pct is not 0:
        decimals = str(min_pct)[::-1].find('.')
        min_pct = min_pct / 100.0
    else:
        decimals = 1
        min_pct = 0.001

    # set up plot
    fig, ax = plt.subplots(1, 1, figsize=figsize, tight_layout=True)
    plt.xlim([1,max_freq])
    plt.ylim([min_pct,1])
    plt.xscale('log')
    plt.yscale('log')
    plt.xlabel('Frequency (Hz)')
    plt.ylabel('Modulation (%)')
    plt.gca().xaxis.set_major_formatter(ScalarFormatter())
    plt.gca().yaxis.set_major_formatter(PercentFormatter(1, decimals=decimals))
    plt.grid(which='both')
    ax.set_axisbelow(True)

    # plot no risk region
    norisk_region = [[1, min_pct], [1, 0.001], [10, 0.001], [100, 0.01], [100, 0.03], [3000, 1], \
        [max_freq, 1], [max_freq, min_pct], ]
    fc_color = 'g'
    if not noriskcolor:
        fc_color = 'gray'
    norisk = plt.Polygon(norisk_region, fc=fc_color, alpha=0.3)
    plt.gca().add_patch(norisk)

    # plot low risk region
    if showyellow:
        lowrisk_region = [[1, 0.001], [1, 0.002], [8, 0.002], [90, 0.025], [90, 0.075], [1200, 1], \
            [3000, 1], [100, 0.03], [100, 0.025], [100, 0.01], [10, 0.001]]
        lowrisk = plt.Polygon(lowrisk_region, fc='y', alpha=0.3)
        plt.gca().add_patch(lowrisk)

    # plot high risk region
    if showred:
        highrisk_region = [[1, 0.002], [8, 0.002], [90, 0.025], [90, 0.075], [1200, 1], [1, 1]]
        highrisk = plt.Polygon(highrisk_region, fc='r', alpha=0.2)
        plt.gca().add_patch(highrisk)

    # plot the data
    markers = itertools.cycle(('o', '^', 's', 'D', 'p', 'P'))
    for pt in data:
        plt.scatter(pt[0], pt[1], label=pt[2], marker=next(markers), alpha=1)
    plt.legend()

    # save the figure if a filename was specified
    if filename:
        plt.savefig(filename, dpi=300)

    # show the plot
    if not suppress:
        plt.show()

def waveform_graph(waveform, figsize:tuple=(8,4), suppress:bool=False, filename:str=None,
                   showstats:bool=True, showstandards:bool=True, num_periods:int=None,
                   fullheight:bool=False, data=None):
    """Plots the time-domain flicker waveform

    Parameters
    ----------
    waveform : Waveform
        The Waveform object
    figsize : tuple
        The (horizontal, vertical) figure size
    suppress : bool
        If True, will not run plt.show()
    filename : str or None
        If specified, will save the file named as such, e.g.: filename='../out/this_graph.png'
    showstats : bool
        If True, will show the flicker frequency, percent, and index on the bottom left of the graph
    showstandards : bool
        If True, will display IEEE, WELL, and JA8 test results on the bottom right of the graph
    num_periods : int or None
        If not None, will truncate the graph to the specified number of periods.
        NOTE: Must be less than the number of periods present in the data
    fullheight : bool
        If True, will set the bottom y limit to zero.
        If False, will display the non-scaled waveform starting from v_min
    data : np.ndarray or None
        If None, will plot the regular waveform
        If ndarray, will plot the array
    """

    # Get the data from the waveform
    if data is None:
        data = waveform.get_data()

    # Create the figure
    fig, ax = plt.subplots(1, 1, figsize=figsize)
    plt.ylabel('Light Output')
    plt.xlabel('Time (ms)')

    # hide the top and right axes
    ax.spines['top'].set_color('none')
    ax.spines['right'].set_color('none')

    # get the number of periods to display
    if num_periods:
        data = waveform.get_n_periods(num_periods=num_periods)

    # scale the x axis to milliseconds
    x_data = data[:,0] * 1000

    # get minimum for y axis scaling
    y_min = waveform.get_v_min() / waveform.get_v_max()

    # display the waveform full height? (xmin=0)
    if fullheight:
        plt.ylim((0,1))
        plt.yticks(np.linspace(0, 1, 6))

        # scale y axis to 0.99 because of strange clipping at 1.0
        y_data = minmax_scale(data[:,1], feature_range=(y_min,0.99))
    else:
        ax.spines['left'].set_smart_bounds(True)

        # scale y axis to (0,1)
        y_data = minmax_scale(data[:,1], feature_range=(y_min,0.99))

    # make the left and bottom axis look cleaner
    ax.spines['bottom'].set_smart_bounds(True)

    # plot
    ax.plot(x_data, y_data)

    # show stats on the graph
    if showstats:
        text(0.02, 0.1, waveform.summary(), ha='left', va='center', transform=ax.transAxes)

    # show standard test results on the graph
    if showstandards:
        std_text = "IEEE 1789: " + waveform.get_ieee_1789_2015() + \
            "\nCalifornia JA8: " +  bool_to_pass_fail(waveform.get_california_ja8_2019()) + \
            "\nWELL v2: " +  bool_to_pass_fail(waveform.get_well_standard_v2())
        text(0.955, 0.1, std_text, ha='right', va='center', transform=ax.transAxes) #, backgroundcolor='silver')

    # save the figure if a filename was specified
    if filename:
        plt.savefig(filename, dpi=300)

    # show the plot
    if not suppress:
        plt.show()

def standards_color(result:str) -> str:
    """Returns colors for decorating standards result labels

    For example, "Pass" returns 'green' and "Fail" returns 'red'

    Parameters
    ----------
    result : str
        The text string to apply color to

    Returns
    -------
    str
        The color, either 'green', 'yellow', or 'red'
    """

    if result is "Pass" or result is "No Risk":
        return 'green'
    elif result is "Low Risk":
        return 'yellow'
    else:
        return 'red'

max_stars_repo_name: yeutterg/beautiful-flicker

max_issues_repo_name: yeutterg/beautiful-flicker

max_forks_repo_name: yeutterg/beautiful-flicker


text:
# Copyright (c) Facebook, Inc. and its affiliates.
#
# This source code is licensed under the MIT license found in the
# LICENSE file in the root directory of this source tree.
import logging
import shutil
import subprocess
import tempfile
from pathlib import Path
from threading import Lock
from typing import Iterable, Optional, Union

import numpy as np
from fasteners import InterProcessLock

from compiler_gym.datasets import Benchmark, BenchmarkSource, Dataset
from compiler_gym.datasets.benchmark import BenchmarkWithSource
from compiler_gym.envs.gcc.gcc import Gcc
from compiler_gym.util.decorators import memoized_property
from compiler_gym.util.runfiles_path import runfiles_path
from compiler_gym.util.shell_format import plural
from compiler_gym.util.truncate import truncate

# The maximum value for the --seed argument to csmith.
UINT_MAX = (2 ** 32) - 1

_CSMITH_BIN = runfiles_path("compiler_gym/third_party/csmith/csmith/bin/csmith")
_CSMITH_INCLUDES = runfiles_path(
    "compiler_gym/third_party/csmith/csmith/include/csmith-2.3.0"
)
_CSMITH_INSTALL_LOCK = Lock()

# TODO(github.com/facebookresearch/CompilerGym/issues/325): This can be merged
# with the LLVM implementation.
class CsmithBenchmark(BenchmarkWithSource):
    """A CSmith benchmark."""

    def __init__(self, *args, **kwargs):
        super().__init__(*args, **kwargs)
        self._src = None

    @classmethod
    def create(cls, uri: str, bitcode: bytes, src: bytes) -> Benchmark:
        """Create a benchmark from paths."""
        benchmark = cls.from_file_contents(uri, bitcode)
        benchmark._src = src  # pylint: disable=protected-access
        return benchmark

    @memoized_property
    def sources(self) -> Iterable[BenchmarkSource]:
        return [
            BenchmarkSource(filename="source.c", contents=self._src),
        ]

    @property
    def source(self) -> str:
        """Return the single source file contents as a string."""
        return self._src.decode("utf-8")

class CsmithDataset(Dataset):
    """A dataset which uses Csmith to generate programs.

    Csmith is a tool that can generate random conformant C99 programs. It is
    described in the publication:

        Yang, Xuejun, Yang Chen, Eric Eide, and John Regehr. "Finding and
        understanding bugs in C compilers." In Proceedings of the 32nd ACM
        SIGPLAN conference on Programming Language Design and Implementation
        (PLDI), pp. 283-294. 2011.

    For up-to-date information about Csmith, see:

        https://embed.cs.utah.edu/csmith/

    Note that Csmith is a tool that is used to find errors in compilers. As
    such, there is a higher likelihood that the benchmark cannot be used for an
    environment and that :meth:`env.reset()
    <compiler_gym.envs.CompilerEnv.reset>` will raise :class:`BenchmarkInitError
    <compiler_gym.datasets.BenchmarkInitError>`.
    """

    def __init__(
        self,
        gcc_bin: Union[Path, str],
        site_data_base: Path,
        sort_order: int = 0,
        csmith_bin: Optional[Path] = None,
        csmith_includes: Optional[Path] = None,
    ):
        """Constructor.

        :param site_data_base: The base path of a directory that will be used to
            store installed files.

        :param sort_order: An optional numeric value that should be used to
            order this dataset relative to others. Lowest value sorts first.

        :param csmith_bin: The path of the Csmith binary to use. If not
            provided, the version of Csmith shipped with CompilerGym is used.

        :param csmith_includes: The path of the Csmith includes directory. If
            not provided, the includes of the Csmith shipped with CompilerGym is
            used.
        """
        super().__init__(
            name="generator://csmith-v0",
            description="Random conformant C99 programs",
            references={
                "Paper": "http://web.cse.ohio-state.edu/~rountev.1/5343/pdf/pldi11.pdf",
                "Homepage": "https://embed.cs.utah.edu/csmith/",
            },
            license="BSD",
            site_data_base=site_data_base,
            sort_order=sort_order,
            benchmark_class=CsmithBenchmark,
        )
        self.gcc_bin = gcc_bin
        self.csmith_bin_path = csmith_bin or _CSMITH_BIN
        self.csmith_includes_path = csmith_includes or _CSMITH_INCLUDES
        self._install_lockfile = self.site_data_path / ".install.LOCK"

    @property
    def size(self) -> int:
        # Actually 2^32 - 1, but practically infinite for all intents and
        # purposes.
        return 0

    @memoized_property
    def gcc(self):
        # Defer instantiation of Gcc from the constructor as it will fail if the
        # given Gcc is not available. Memoize the result as initialization is
        # expensive.
        return Gcc(bin=self.gcc_bin)

    def benchmark_uris(self) -> Iterable[str]:
        return (f"{self.name}/{i}" for i in range(UINT_MAX))

    def benchmark(self, uri: str) -> CsmithBenchmark:
        return self.benchmark_from_seed(int(uri.split("/")[-1]))

    def _random_benchmark(self, random_state: np.random.Generator) -> Benchmark:
        seed = random_state.integers(UINT_MAX)
        return self.benchmark_from_seed(seed)

    @property
    def installed(self) -> bool:
        return super().installed and (self.site_data_path / "includes").is_dir()

    def install(self) -> None:
        super().install()

        if self.installed:
            return

        with _CSMITH_INSTALL_LOCK, InterProcessLock(self._install_lockfile):
            if (self.site_data_path / "includes").is_dir():
                return

            # Copy the Csmith headers into the dataset's site directory path because
            # in bazel builds this includes directory is a symlink, and we need
            # actual files that we can use in a docker volume.
            shutil.copytree(
                self.csmith_includes_path,
                self.site_data_path / "includes.tmp",
            )
            # Atomic directory rename to prevent race on install().
            (self.site_data_path / "includes.tmp").rename(
                self.site_data_path / "includes"
            )

    def benchmark_from_seed(
        self, seed: int, max_retries: int = 3, retry_count: int = 0
    ) -> CsmithBenchmark:
        """Get a benchmark from a uint32 seed.

        :param seed: A number in the range 0 <= n < 2^32.

        :return: A benchmark instance.

        :raises OSError: If Csmith fails.

        :raises BenchmarkInitError: If the C program generated by Csmith cannot
            be lowered to LLVM-IR.
        """
        if retry_count >= max_retries:
            raise OSError(
                f"Csmith failed after {retry_count} {plural(retry_count, 'attempt', 'attempts')} "
                f"with seed {seed}"
            )

        self.install()

        # Run csmith with the given seed and pipe the output to clang to
        # assemble a bitcode.
        self.logger.debug("Exec csmith --seed %d", seed)
        csmith = subprocess.Popen(
            [str(self.csmith_bin_path), "--seed", str(seed)],
            stdout=subprocess.PIPE,
            stderr=subprocess.PIPE,
        )

        # Generate the C source.
        src, stderr = csmith.communicate(timeout=300)
        if csmith.returncode:
            try:
                stderr = "\n".join(
                    truncate(stderr.decode("utf-8"), max_line_len=200, max_lines=20)
                )
                logging.warning("Csmith failed with seed %d: %s", seed, stderr)
            except UnicodeDecodeError:
                # Failed to interpret the stderr output, generate a generic
                # error message.
                logging.warning("Csmith failed with seed %d", seed)
            return self.benchmark_from_seed(
                seed, max_retries=max_retries, retry_count=retry_count + 1
            )

        # Pre-process the source.
        with tempfile.TemporaryDirectory() as tmpdir:
            src_file = f"{tmpdir}/src.c"
            with open(src_file, "wb") as f:
                f.write(src)

            preprocessed_src = self.gcc(
                "-E",
                "-I",
                str(self.site_data_path / "includes"),
                "-o",
                "-",
                src_file,
                cwd=tmpdir,
                timeout=60,
                volumes={
                    str(self.site_data_path / "includes"): {
                        "bind": str(self.site_data_path / "includes"),
                        "mode": "ro",
                    }
                },
            )

        return self.benchmark_class.create(
            f"{self.name}/{seed}", preprocessed_src.encode("utf-8"), src
        )

max_stars_repo_name: AkillesAILimited/CompilerGym

max_issues_repo_name: AkillesAILimited/CompilerGym

max_forks_repo_name: AkillesAILimited/CompilerGym


text:
import numpy as np
import random
import os
import sys

sys.path.append("../")
import pyedflib

# from constants import INCLUDED_CHANNELS, FREQUENCY
from scipy.fftpack import fft
from scipy.signal import resample

INCLUDED_CHANNELS = [
    "EEG FP1",
    "EEG FP2",
    "EEG F3",
    "EEG F4",
    "EEG C3",
    "EEG C4",
    "EEG P3",
    "EEG P4",
    "EEG O1",
    "EEG O2",
    "EEG F7",
    "EEG F8",
    "EEG T3",
    "EEG T4",
    "EEG T5",
    "EEG T6",
    "EEG FZ",
    "EEG CZ",
    "EEG PZ",
]

FREQUENCY = 200

def computeFFT(signals, n):
    """
    Args:
        signals: EEG signals, (number of channels, number of data points)
        n: length of positive frequency terms of fourier transform
    Returns:
        FT: log amplitude of FFT of signals, (number of channels, number of data points)
        P: phase spectrum of FFT of signals, (number of channels, number of data points)
    """
    # fourier transform
    fourier_signal = fft(signals, n=n, axis=-1)  # FFT on the last dimension

    # only take the positive freq part
    idx_pos = int(np.floor(n / 2))
    fourier_signal = fourier_signal[:, :idx_pos]
    amp = np.abs(fourier_signal)
    amp[amp == 0.0] = 1e-8  # avoid log of 0

    FT = np.log(amp)
    P = np.angle(fourier_signal)

    return FT, P

def get_swap_pairs(channels):
    """
    Swap select adjacenet channels
    Args:
        channels: list of channel names
    Returns:
        list of tuples, each a pair of channel indices being swapped
    """
    swap_pairs = []
    if ("EEG FP1" in channels) and ("EEG FP2" in channels):
        swap_pairs.append((channels.index("EEG FP1"), channels.index("EEG FP2")))
    if ("EEG Fp1" in channels) and ("EEG Fp2" in channels):
        swap_pairs.append((channels.index("EEG Fp1"), channels.index("EEG Fp2")))
    if ("EEG F3" in channels) and ("EEG F4" in channels):
        swap_pairs.append((channels.index("EEG F3"), channels.index("EEG F4")))
    if ("EEG F7" in channels) and ("EEG F8" in channels):
        swap_pairs.append((channels.index("EEG F7"), channels.index("EEG F8")))
    if ("EEG C3" in channels) and ("EEG C4" in channels):
        swap_pairs.append((channels.index("EEG C3"), channels.index("EEG C4")))
    if ("EEG T3" in channels) and ("EEG T4" in channels):
        swap_pairs.append((channels.index("EEG T3"), channels.index("EEG T4")))
    if ("EEG T5" in channels) and ("EEG T6" in channels):
        swap_pairs.append((channels.index("EEG T5"), channels.index("EEG T6")))
    if ("EEG O1" in channels) and ("EEG O2" in channels):
        swap_pairs.append((channels.index("EEG O1"), channels.index("EEG O2")))

    return swap_pairs

def random_augmentation(EEG_seq):
    """
    Random augmentation of EEG sequence by randomly swapping channels
    and randomly scaling the signals.
    Args:
        EEG_seq: shape (seq_length, num_nodes, num_data_point)
    Returns:
        augmented signals, same shape as input
    """
    for pair in get_swap_pairs(INCLUDED_CHANNELS):
        if random.choice([True, False]):
            EEG_seq[:, [pair[0], pair[1]], :] = EEG_seq[:, [pair[1], pair[0]], :]
    if random.choice([True, False]):
        EEG_seq = EEG_seq * np.random.uniform(0.8, 1.2)
    return EEG_seq

def getOrderedChannels(file_name, verbose, labels_object, channel_names):
    labels = list(labels_object)
    for i in range(len(labels)):
        labels[i] = labels[i].split("-")[0]

    ordered_channels = []
    for ch in channel_names:
        try:
            ordered_channels.append(labels.index(ch))
        except:
            if verbose:
                print(file_name + " failed to get channel " + ch)
            raise Exception("channel not match")
    return ordered_channels

def getChannelIndices(labelsObject, channel_names):
    """
    labelsObject: labels object from eeghdf file
    channel_names: list of channel names
    """
    labels = labelsObject
    channelIndices = [labels.index(ch) for ch in channel_names]
    return channelIndices

def getSeizureTimes(file_name):
    """
    Args:
        file_name: file name of .edf file etc.
    Returns:
        seizure_times: list of times of seizure onset in seconds
    """
    tse_file = file_name.split(".edf")[0] + ".tse_bi"

    seizure_times = []
    with open(tse_file) as f:
        for line in f.readlines():
            if "seiz" in line:  # if seizure
                # seizure start and end time
                seizure_times.append(
                    [
                        float(line.strip().split(" ")[0]),
                        float(line.strip().split(" ")[1]),
                    ]
                )
    return seizure_times

def getEDFsignals(edf):
    n = edf.signals_in_file
    samples = edf.getNSamples()[0]
    signals = np.zeros((n, samples))
    for i in range(n):
        try:
            signals[i, :] = edf.readSignal(i)
        except:
            pass
    return signals

def resampleData(signals, to_freq=200, window_size=4):
    """
    Resample signals from its original sampling freq to another freq
    Args:
        signals: EEG signal slice, (num_channels, num_data_points)
        to_freq: Re-sampled frequency in Hz
        window_size: time window in seconds
    Returns:
        resampled: (num_channels, resampled_data_points)
    """
    num = int(to_freq * window_size)
    resampled = resample(signals, num=num, axis=1)
    return resampled

class StandardScaler:
    """
    Standardize the input
    """

    def __init__(self, mean, std):
        self.mean = mean  # (1,num_nodes,1)
        self.std = std  # (1,num_nodes,1)

    def transform(self, data):
        return (data - self.mean) / self.std

    def inverse_transform(self, data):
        """
        Masked inverse transform
        Args:
            data: data for inverse scaling
            is_tensor: whether data is a tensor
            device: device
        """
        return data * self.std + self.mean

max_stars_repo_name: dumpmemory/state-spaces

max_issues_repo_name: dumpmemory/state-spaces

max_forks_repo_name: MikeOwino/state-spaces

reason: import numpy,from scipy


text:
from __future__ import division
import numpy as np
from DSGC_SAC import *

# Parameters ##################################################################
n_angles = 4  # Number of angles to plot
nsyn = 8  # Synaptic connection to plot (see synM)
tstop = 1000  # Simulation duration
SAC_excGmax = 4.25e-5  # Stim amplitude (G_max excitatory) for SACs
bar_speed = 1.61
excGmax = 1 * SAC_excGmax  # Stim amplitude for DSGC; note factor
gabaGmax = 0.46 * s.gabaGmax  # Inhibitory weight w/r to SAC
olap_thresh = 1000  # Threshold of overlap area for SAC/DSGC connection
ofactor = 1  # Weight of overlap area for synapse
###############################################################################

diag = np.sqrt(dsgc_x ** 2 + dsgc_y ** 2)
if n_angles == 12:
    angles = (0, 30, 60, 90, 120, 150, 180, 210, 240, 270, 300, 330)
    d_inits = (0, 0, 0, 0, 2 * dsgc_x, 2 * dsgc_x, 2 * dsgc_x, 2 * diag,
               2 * diag, 2 * dsgc_y, 2 * dsgc_y, 2 * dsgc_y)
    plot_indexs = (6, 3, 3, 2, 1, 1, 4, 7, 7, 8, 9, 9)
elif n_angles == 8:
    angles = (0, 45, 90, 135, 180, 225, 270, 315)
    d_inits = (0, 0, 0, diag, 2 * dsgc_x, 2 * diag, 2 * dsgc_y, diag)
    plot_indexs = (6, 3, 2, 1, 4, 7, 8, 9)
elif n_angles == 4:
    angles = (0, 180, 90, 270)
    d_inits = (0, 2 * dsgc_x, 0, 2 * dsgc_y)
    plot_indexs = (6, 4, 2, 8)
elif n_angles == 2:
    angles = (180, 0)
    d_inits = (2 * dsgc_x, 0)
    plot_indexs = (4, 6)
elif n_angles == 1:
    angles = (180, )
    d_inits = (2 * dsgc_x, )
    plot_indexz = (4, )

nspk = []
fig_v = plt.figure()
fig_s1 = plt.figure()
fig_s2 = plt.figure()
fig_s3 = plt.figure()

gBPsyn = bipolar_DSGC_synapses(s)
GABAsynpos, synM = SAC_DSGC_synapses(s,
                                     thresh=olap_thresh,
                                     overlap_factor=ofactor)

for (ang, d_init, plti) in zip(angles[::-1], d_inits[::-1], plot_indexs[::-1]):
    stim_param['bar_x_init'] = d_init
    stim_param['bar_angle'] = ang
    stim_param['bar_speed'] = bar_speed
    # Setting DSGC inputs #####################################################
    # This is for DSGC inputs. Note that SAC syanpses are set and defined in
    # SAC_network
    dscell.setv(-75.5, -75.5, -79)
    gBPsyn = s.topo.excGset(gBPsyn, gmax=excGmax, gmin=excGmin,
                            synapse_type=synapse_type, tau=tau)
    s.stim.set_stimulus(gBPsyn, stimulus_type, t_inf, synapse_type,
                        **stim_param)
    s.topo.gabaGset(gSACsyn, synM, gabaGmax, gabaGmin, k1, k2, th1, th2, e_rev)

    # Recording vectors #######################################################
    timevec = s.h.Vector()
    timevec.from_python(np.arange(0, tstop, sampinvl))
    timerec = s.h.Vector()
    timerec.record(s.h._ref_t, timevec)
    dsgcV = [s.h.Vector(tstop / sampinvl) for i in range(7)]
    dsgcV[0].record(dscell.dend1[0](0.5)._ref_v, timevec)
    dsgcV[1].record(dscell.dend1[2](0.5)._ref_v, timevec)
    dsgcV[2].record(dscell.dend1[4](0.5)._ref_v, timevec)
    dsgcV[3].record(dscell.dend1[6](0.5)._ref_v, timevec)
    dsgcV[4].record(dscell.soma(0.75)._ref_v, timevec)
    dsgcV[5].record(dscell.axon(0.07)._ref_v, timevec)
    dsgcV[6].record(dscell.axon(0.99)._ref_v, timevec)

    inh_synI = [s.h.Vector(tstop / sampinvl) for i in range(len(gSACsyn))]
    for (r, syn) in zip(inh_synI, gSACsyn):
        r.record(syn._ref_i, timevec)

    exc_synI = [s.h.Vector(tstop / sampinvl) for i in range(len(gBPsyn))]
    for (r, syn) in zip(exc_synI, gBPsyn):
        r.record(syn[0]._ref_i, timevec)

    vpre = [s.h.Vector(tstop / sampinvl) for i in range(len(gSACsyn))]
    for (r, syn) in zip(vpre, gSACsyn):
        r.record(syn._ref_vpre, timevec)

    sc2 = [s.h.Vector(tstop / sampinvl) for i in range(len(gSACsyn))]
    for (r, syn) in zip(sc2, gSACsyn):
        r.record(syn._ref_sc2, timevec)

    apcount = s.h.APCount(dscell.axon(.9))

    # ginh = [s.h.Vector(tstop / sampinvl) for i in range(len(gSACsyn))]
    # for (r, syn) in zip(ginh, gSACsyn):
    #     r.record(syn._ref_g, timevec)

    # iclamp = insert_soma_stim()
    # pt.bipolar_SAC_map(s, '', dscell, BPpos, GABAsynpos)
    res, amac_vecs, synapse_vecs = s.main(excmax=SAC_excGmax,
                                          gabaI=1, main=1,
                                          stim_param=stim_param,
                                          tstop=tstop, v_init=SAC_v_init)
    nspk.append(apcount.n)

    plt.figure(fig_v.number)
    plt.subplot(3, 3, plti)
    pltleg = 8
    for v in dsgcV:
        plt.plot(timerec, v)
    if plti == pltleg:
        if dsgc_type in ('DSGC', 'DSGC2'):
            plt.legend(('dend1[0]', 'dend1[2]', 'dend1[4]',
                        'dend1[6]', 'soma', 'AIS', 'axon'), loc='best')
        elif dsgc_type == 'DSGC_Poleg':
            plt.legend(('soma', 'dend[0]', 'dend[20]', 'dend[40]',
                        'dend[60]'), loc=7)
        plt.xlabel('Time (ms)')
        plt.ylabel('Transmembrane Voltage (mV)')
    plt.ylim((-80, -30))
    # plt.subplot(212)
    # for a in dsgcA:
    #    plt.plot(timerec,a)
    # plt.legend(('m_soma','h_soma','m_AIS','h_AIS'))
    plt.figure(fig_s1.number)
    plt.subplot(3, 3, plti)
    if plti == pltleg:
        plt.title('GABA Synapse #' + str(nsyn))
        plt.ylabel('Preterminal Voltage (mV)')
        plt.xlabel('Time (ms)')
    plt.plot(timevec, vpre[nsyn])
    plt.ylim((-65, -40))

    for ns in range(0, len(gSACsyn)):
        plt.figure(fig_s2.number)
        plt.subplot(3, 3, plti)
        g = gSACsyn[ns].minG + (
            gSACsyn[ns].maxG - gSACsyn[ns].minG) * sc2[ns].as_numpy()
        i2 = g * (dsgcV[0].as_numpy() - e_rev)
        plt.plot(timevec, 1000 * i2, timevec, 1000 * exc_synI[nsyn].as_numpy())
        if plti == pltleg:
            plt.title('GABA Synapse #' + str(nsyn))
            plt.ylabel('Synaptic Current ($\mu$A)')
            plt.xlabel('Time (ms)')
        # plt.ylim((0, .20))

        plt.figure(fig_s3.number)
        plt.subplot(3, 3, plti)
        plt.plot(timevec, g * 1000)
        if plti == pltleg:
            plt.title('GABA synapse #' + str(nsyn))
            plt.ylabel('Synaptic Conductance (mS)')
            plt.xlabel('Time (ms)')
        # plt.ylim((0, .150))

dsi = s.tl.DSI(angles, nspk)


text:
import os
import random
import numpy as np
from skimage import io, transform
from PIL import Image

Image.MAX_IMAGE_PIXELS = None

def make_dir(path):
    if not os.path.exists(path):
        os.makedirs(path)

root_dir = r'C:\Users\bunny\Desktop\data'
source_img_path = r"C:\Users\bunny\Desktop\merge.tif"
source_img = io.imread(source_img_path)
source_mask_path = r"C:\Users\bunny\Desktop\merge_db_mask_.tif"
source_mask = io.imread(source_mask_path)
split_func = ['train', 'test']
split_ratio = 0.8
img_count = 5000
crop_size = 512
functions = ['image', 'mask']

sp_line = int(source_img.shape[0] * split_ratio)
sp_info = {
    'train': {'img': {'image': source_img[:sp_line],
                      'mask': source_mask[:sp_line], },
              'count': int(img_count * split_ratio), },
    'test': {'img': {'image': source_img[sp_line:],
                     'mask': source_mask[sp_line:], },
             'count': int(img_count * (1 - split_ratio))},
}

for sp in sp_info:
    print(sp)
    for function in sp_info[sp]['img']:
        output_dir = os.path.join(root_dir, sp, function)
        make_dir(output_dir)
    image_count = 0

    x_max, y_max = sp_info[sp]['img']['image'].shape[:2]
    while image_count < sp_info[sp]['count']:
        startpoint = [random.randint(0, x_max - crop_size),
                      random.randint(0, y_max - crop_size)]
        end_point = [point + crop_size for point in startpoint]
        for function in sp_info[sp]['img']:
            output_img = sp_info[sp]['img'][function][startpoint[0]: end_point[0], startpoint[1]: end_point[1]]
            loc = '_'.join([str(point) for point in startpoint + end_point])
            output_path = os.path.join(root_dir, sp, function, loc + '.png')
            io.imsave(output_path, output_img, check_contrast=False)
        image_count += 1
        print('\rImage done: {}'.format(image_count), end='')
    print()

max_stars_repo_name: juyingnan/siie

max_issues_repo_name: juyingnan/siie

max_forks_repo_name: juyingnan/siie


text:
'''This script demonstrates how to build a variational autoencoder
with Keras and deconvolution layers.

Reference: "Auto-Encoding Variational Bayes" https://arxiv.org/abs/1312.6114

Using Tensorflow Dataset API and Horovod.

original implementation:
https://github.com/fchollet/keras/blob/master/examples/variational_autoencoder_deconv.py

run:

    # 2 GPUs example
    mpirun -np 2 --map-by ppr:4:socket python \
        variational_autoencoder_deconv_tfqueue_horovod.py --epochs=30

    # 4 GPUs example
    TMPDIR=/tmp mpirun --report-bindings --bind-to none --map-by slot -np 4 \
      python ./examples/variational_autoencoder/variational_autoencoder_deconv_tfdataset_horovod.py \
      --epochs=30  --nranks_per_gpu=1

'''  # noqa
from __future__ import print_function
import sys
import argparse as ap

import numpy as np

# try:
#     import Tkinter  # @UnusedImport
#     import matplotlib.pyplot as plt
# except ImportError:
#     import matplotlib
#     matplotlib.use('Agg')
#     import matplotlib.pyplot as plt
import matplotlib
try:
    matplotlib.use('Agg')
except Exception:
    raise
import matplotlib.pyplot as plt

from scipy.stats import norm

import tensorflow as tf
try:
    Dataset = tf.data.Dataset
except Exception:
    from tensorflow.contrib.data import Dataset

import horovod.tensorflow as hvd
# import horovod.keras as hvd_keras

try:
    # Initialize Horovod.
    hvd.init()
except Exception:
    raise

from keras.layers import Input
from keras.models import Model
from keras import backend as K
from keras.datasets import mnist
from keras.optimizers import TFOptimizer
# from keras.optimizers import RMSprop

from keras_exp.callbacks.timing import BatchTiming, SamplesPerSec

from vae_common import (
    CustomFormatter, make_shared_layers_dict,
    make_vae, get_encoded, get_decoded)

class SpeedupOpts(object):
    epoch = 'epoch'
    imgspersec = 'imgspersec'

def parser_(desc):
    parser = ap.ArgumentParser(description=desc,
                               formatter_class=CustomFormatter)

    parser.add_argument(
        '--epochs', type=int, default=5,
        help='Number of epochs to run training for. Default: %(default)s')

    parser.add_argument(
        '--nranks_per_gpu', type=int, default=1,
        help='S|Number of ranks to run on each GPUs. Use this parameter to\n'
        'oversubscribe a GPU. When oversubscribing a GPU use in combination\n'
        'with MPS (multi-process service). Default: %(default)s')

    parser.add_argument(
        '--speedup', action='store', nargs='?', type=str.lower,
        const=SpeedupOpts.epoch, default=SpeedupOpts.epoch,
        choices=[SpeedupOpts.epoch, SpeedupOpts.imgspersec],
        help='S|Speedup the fit to run each epoch faster or overall images\n'
        'per seconds faster. With the "imgspersec" option each rank\n'
        'fits the fullset and with epoch option only a subset. Loss\n'
        'convergence stability might be better with imgspersec option.\n'
        'Default: %(default)s')

    args = parser.parse_args()

    return args

def main(argv=None):
    '''
    '''
    main.__doc__ = __doc__
    argv = sys.argv if argv is None else sys.argv.extend(argv)
    desc = main.__doc__  # .format(os.path.basename(__file__))
    # CLI parser
    args = parser_(desc)

    nranks_per_gpu = args.nranks_per_gpu
    local_rank = hvd.local_rank()
    gpu_local_rank = local_rank // nranks_per_gpu
    print('local_rank, GPU_LOCAL_RANK: {}, {}'.format(
        local_rank, gpu_local_rank))

    # Pin GPU to be used to process local rank (one GPU per process)
    config = tf.ConfigProto()
    config.gpu_options.allow_growth = True
    # config.gpu_options.visible_device_list = str(hvd.local_rank())
    config.gpu_options.visible_device_list = str(gpu_local_rank)
    K.set_session(tf.Session(config=config))

    # input image dimensions
    img_rows, img_cols, img_chns = 28, 28, 1
    # number of convolutional filters to use
    filters = 64
    # convolution kernel size
    num_conv = 3

    hvdsize = hvd.size()

    batch_size = 128  # 100
    if K.image_data_format() == 'channels_first':
        original_img_size = (img_chns, img_rows, img_cols)
    else:
        original_img_size = (img_rows, img_cols, img_chns)
    latent_dim = 2
    intermediate_dim = 128
    epsilon_std = 1.0
    epochs = args.epochs  # 5

    # train the VAE on MNIST digits
    (x_train, _), (x_test, y_test) = mnist.load_data()

    x_train = x_train.astype('float32') / 255.
    x_train = x_train.reshape((x_train.shape[0],) + original_img_size)
    x_test = x_test.astype('float32') / 255.
    x_test = x_test.reshape((x_test.shape[0],) + original_img_size)

    if hvd.rank() == 0:
        print('x_train.shape:', x_train.shape)

    train_samples = x_train.shape[0]
    # steps_per_epoch = train_samples // batch_size // hvdsize
    speedupopt = args.speedup
    if speedupopt == SpeedupOpts.imgspersec:
        steps_per_epoch = train_samples // batch_size
    else:
        steps_per_epoch = int(round(
            float(train_samples) / batch_size / hvdsize + 0.5))

    # Create the dataset and its associated one-shot iterator.
    buffer_size = 10000
    dataset = Dataset.from_tensor_slices(x_train)
    dataset = dataset.repeat()
    dataset = dataset.shuffle(buffer_size)
    dataset = dataset.batch(batch_size)
    iterator = dataset.make_one_shot_iterator()
    x_train_batch = iterator.get_next()

    ldict = make_shared_layers_dict(
        img_chns, img_rows, img_cols, batch_size, filters,
        num_conv, intermediate_dim, latent_dim, epsilon_std)
    # ldict is a dictionary that holds all layers. Since these layers are
    # instantiated once, they are shared amongs vae, encoder, and generator.

    x = Input(tensor=x_train_batch)
    vae = make_vae(ldict, x)
    # :  :type vae: Model

    lr = 0.001  # * hvdsize
    opt = tf.train.RMSPropOptimizer(lr)
    # Add Horovod Distributed Optimizer.
    opt = hvd.DistributedOptimizer(opt)  # , use_locking=True)
    opt = TFOptimizer(opt)

    # opt = RMSprop(lr)
    # Add Horovod Distributed Optimizer.
    # opt = hvd_keras.DistributedOptimizer(opt)  # , use_locking=True)

    vae.compile(optimizer=opt, loss=None)
    if hvd.rank() == 0:
        vae.summary()

    callbacks = []
    if hvd.rank() == 0:
        callbacks += [BatchTiming(), SamplesPerSec(batch_size * hvdsize)]

    sess = K.get_session()
    sess.run(hvd.broadcast_global_variables(0))

    # Fit the model using data from the TF data tensors.
    vae.fit(steps_per_epoch=steps_per_epoch, epochs=epochs,
            callbacks=callbacks)

    if hvd.rank() == 0:
        x = Input(shape=original_img_size)
        vae_val = make_vae(ldict, x)
        vae_val.compile(optimizer=opt, loss=None)
        loss = vae_val.evaluate(x=x_test, y=None, batch_size=batch_size)
        print('\n\nVAE VALIDATION LOSS: {}'.format(loss))

        x = Input(shape=original_img_size)
        z_mean, _ = get_encoded(ldict, x)
        encoder = Model(x, z_mean)
        # :  :type encoder: Model

        decoder_input = Input(shape=(latent_dim,))
        x_decoded_mean_squash = get_decoded(ldict, decoder_input)
        generator = Model(decoder_input, x_decoded_mean_squash)
        # :  :type generator: Model

        # display a 2D plot of the digit classes in the latent space
        x_test_encoded = encoder.predict(x_test, batch_size=batch_size)
        plt.figure(figsize=(6, 6))
        plt.scatter(x_test_encoded[:, 0], x_test_encoded[:, 1], c=y_test)
        plt.colorbar()
        # plt.show()
        plt.savefig('vae_scatter.ps')
        plt.close()

        # display a 2D manifold of the digits
        n = 15  # figure with 15x15 digits
        digit_size = 28
        figure = np.zeros((digit_size * n, digit_size * n))
        # Linearly spaced coordinates on the unit square were transformed
        # through the inverse CDF (ppf) of the Gaussian
        # To produce values of the latent variables z, since the prior of the
        # latent space is Gaussian
        grid_x = norm.ppf(np.linspace(0.05, 0.95, n))
        grid_y = norm.ppf(np.linspace(0.05, 0.95, n))

        for i, yi in enumerate(grid_x):
            for j, xi in enumerate(grid_y):
                z_sample = np.array([[xi, yi]])
                z_sample = np.tile(z_sample, batch_size).reshape(batch_size, 2)
                x_decoded = generator.predict(z_sample, batch_size=batch_size)
                digit = x_decoded[0].reshape(digit_size, digit_size)
                figure[i * digit_size: (i + 1) * digit_size,
                       j * digit_size: (j + 1) * digit_size] = digit

        plt.figure(figsize=(10, 10))
        plt.imshow(figure, cmap='Greys_r')
        # plt.show()
        plt.savefig('vae_digit.ps')
        plt.close()

    K.clear_session()

if __name__ == '__main__':
    main()


text:
import torch
import torch.nn as nn
import torch.nn.functional as F
from torch.autograd import Function
from torch.nn.modules.module import Module
from torch.nn import init
from collections import OrderedDict
import numpy as np
import correlation_cuda
from torch.nn import Parameter
from torch.optim import lr_scheduler
from lib.warp import backWarp
from lib.visualTool import viz

# for Base option -------------------------------------------------------------------------------

class Interp(nn.Module):
    def __init__(self, scale=None, size=None):
        super(Interp, self).__init__()
        self.scale = scale
        self.size = size

    def forward(self, x):
        return F.interpolate(x, size=self.size, scale_factor=self.scale, mode='bilinear', align_corners=True)

class normalize(nn.Module):
    def __init__(self, dim=1):
        super(normalize, self).__init__()
        self.norm = F.normalize
        self.dim = dim

    def forward(self, tensor: torch.Tensor):
        output = self.norm(tensor, dim=self.dim)
        return output

class Unfold(nn.Module):
    def __init__(self, scale=1, kSzie=None, padding=None):
        super(Unfold, self).__init__()
        self.shape = None
        self.stride = scale
        self.kSzie = (kSzie, kSzie) if kSzie is not None else (3 * self.stride, 3 * self.stride)
        self.padding = padding if padding is not None else self.stride

    def forward(self, inTensor, rev=False):
        if not rev:
            self.shape = inTensor.shape
            outTensor = F.unfold(inTensor,
                                 kernel_size=self.kSzie,
                                 padding=self.padding,
                                 stride=self.stride).reshape(self.shape[0],
                                                             -1,
                                                             self.shape[-2] // self.stride,
                                                             self.shape[-1] // self.stride)
        else:
            inTensor = inTensor.reshape(self.shape[0], -1, self.shape[-2] * self.shape[-1] // self.stride ** 2)
            outTensor = F.fold(inTensor,
                               output_size=(self.shape[-2], self.shape[-1]),
                               kernel_size=self.kSzie,
                               padding=self.padding, stride=self.stride) / (self.kSzie[0] * self.kSzie[1])
        return outTensor

class Swish(nn.Module):
    def __init__(self):
        super(Swish, self).__init__()

    def forward(self, x):
        return x * F.sigmoid(x)

def getPad(in_, kernel_size, stride, dilation=1):
    out_ = np.ceil(float(in_) / stride)
    return int(((out_ - 1) * stride + dilation * (kernel_size - 1) + 1 - in_) / 2)

def getNorm(name: str, num_features: int = 0, affine=True):
    if name.lower() == 'instance':
        return nn.InstanceNorm2d(num_features=num_features, affine=affine)
    if name.lower() == 'lrn':
        return nn.LocalResponseNorm(8)
    if name.lower() == 'bn':
        return nn.BatchNorm2d(num_features=num_features, affine=affine)
    if name.lower() == 'identity':
        return nn.Identity()
    if name.lower() == 'group':
        return torch.nn.GroupNorm(num_groups=4, num_channels=num_features, affine=affine)

def getAct(name: str, num_parameters: int = 0, inplace=True):
    if name.lower() == 'relu':
        return nn.ReLU(inplace=inplace)
    if name.lower() == 'prelu':
        return nn.PReLU(num_parameters)
    if name.lower() == 'leakyrelu':
        return nn.LeakyReLU(negative_slope=0.1, inplace=inplace)
    if name.lower() == 'swish':
        return Swish()

class ConvBlock(nn.Module):
    def __init__(self, in_channels, out_channels, kernel_size=3,
                 stride=1, dilation: tuple = (1, 1), groups=1, padModel='replicate', cfg=None):
        super(ConvBlock, self).__init__()

        if isinstance(kernel_size, tuple):
            kernel_size = kernel_size
        else:
            kernel_size = (kernel_size, kernel_size)
        self.bias = False if 'identity' == cfg.netNorm else True
        self.padModel = padModel

        self.padList = [(kernel_size[1] - 1) // 2, (kernel_size[1] - 1) // 2,
                        (kernel_size[0] - 1) // 2, (kernel_size[0] - 1) // 2]

        self.conv = nn.Sequential(nn.ReplicationPad2d(self.padList),
                                  nn.Conv2d(in_channels=in_channels, out_channels=out_channels, kernel_size=kernel_size,
                                            stride=stride, dilation=dilation, groups=groups, bias=self.bias)
                                  )
        self.norm2d = getNorm(cfg.netNorm, num_features=out_channels)
        self.activate = getAct(cfg.netActivate, num_parameters=out_channels)

    def forward(self, x):

        # x = F.pad(input=x, pad=self.padList, mode=self.padModel)
        x = self.conv(x)
        x = self.norm2d(x)
        x = self.activate(x)

        return x

# end BaseOption -----------------------------------------------------------------------------------

# for pwcNet ---------------------------------------------------------------------------------------
class CorrelationFunction(Function):

    @staticmethod
    def forward(ctx, *agrs):
        input1 = agrs[0]
        input2 = agrs[1]

        ctx.pad_size = agrs[2]

        ctx.kernel_size = agrs[3]
        ctx.max_displacement = agrs[4]
        ctx.stride1 = agrs[5]
        ctx.stride2 = agrs[6]
        ctx.corr_multiply = agrs[7]

        ctx.save_for_backward(input1, input2)

        with torch.cuda.device_of(input1):
            rbot1 = input1.new()
            rbot2 = input2.new()
            output = input1.new()

            correlation_cuda.forward(
                input1,
                input2,
                rbot1,
                rbot2,
                output,
                ctx.pad_size,
                ctx.kernel_size,
                ctx.max_displacement,
                ctx.stride1,
                ctx.stride2,
                ctx.corr_multiply)

        return output

    @staticmethod
    def backward(ctx, *grad_output):
        grad_output = grad_output[0]
        input1, input2 = ctx.saved_tensors

        with torch.cuda.device_of(input1):
            rbot1 = input1.new()
            rbot2 = input2.new()

            grad_input1 = input1.new()
            grad_input2 = input2.new()

            correlation_cuda.backward(
                input1,
                input2,
                rbot1,
                rbot2,
                grad_output,
                grad_input1,
                grad_input2,
                ctx.pad_size,
                ctx.kernel_size,
                ctx.max_displacement,
                ctx.stride1,
                ctx.stride2,
                ctx.corr_multiply)

        return grad_input1, grad_input2, None, None, None, None, None, None

class Correlation(Module):
    def __init__(
            self,
            pad_size=0,
            kernel_size=0,
            max_displacement=0,
            stride1=1,
            stride2=2,
            corr_multiply=1):
        super(Correlation, self).__init__()
        self.pad_size = pad_size
        self.kernel_size = kernel_size
        self.max_displacement = max_displacement
        self.stride1 = stride1
        self.stride2 = stride2
        self.corr_multiply = corr_multiply

    def forward(self, input1, input2):
        Correlation = CorrelationFunction.apply
        result = Correlation(
            input1,
            input2,
            self.pad_size,
            self.kernel_size,
            self.max_displacement,
            self.stride1,
            self.stride2,
            self.corr_multiply)

        return result

def conv(in_planes, out_planes, kernel_size=3, stride=1, padding=1, dilation=1):
    return nn.Sequential(
        nn.Conv2d(in_planes, out_planes, kernel_size=kernel_size, stride=stride,
                  padding=padding, dilation=dilation, bias=True),
        nn.LeakyReLU(0.1, inplace=True))

# end pwcNet --------------------------------------------------------------------------------------------------

# for fuse Net----------------------------------------------------------------------------------------------
class BaseNet(nn.Module):
    """
    init pretrained weights with weight file
    init other weights with a certain random pattern
    """

    def __init__(self, netInitType: str = 'xavier', netInitGain=0.2, netNormAffine=True):
        super(BaseNet, self).__init__()
        self.netInitType = netInitType
        self.netInitGain = netInitGain
        self.netNormAffine = netNormAffine

    def forward(self, *input):
        raise NotImplementedError

    def randomInitNet(self):
        # init all weights by pre-defined pattern firstly
        for m in self.modules():
            if any([isinstance(m, nn.Conv2d), isinstance(m, nn.ConvTranspose2d), isinstance(m, nn.Linear)]):
                if self.netInitType == 'normal':
                    init.normal_(m.weight, 0.0, self.netInitGain)
                elif self.netInitType == 'xavier':
                    init.xavier_normal_(m.weight, gain=self.netInitGain)
                elif self.netInitType == 'kaiming':
                    init.kaiming_normal_(m.weight, a=0.2, mode='fan_in', nonlinearity='leaky_relu')
                elif self.netInitType == 'orthogonal':
                    init.orthogonal_(m.weight, gain=self.netInitGain)
                elif self.netInitType == 'default':
                    pass
                if m.bias is not None:
                    init.constant_(m.bias, 0.0)
            elif any([isinstance(m, nn.InstanceNorm2d), isinstance(m, nn.LocalResponseNorm),
                      isinstance(m, nn.BatchNorm2d), isinstance(m, nn.GroupNorm)]) and self.netNormAffine:
                try:
                    init.constant_(m.weight, 1.0)
                    init.constant_(m.bias, 0.0)
                except Exception as e:
                    pass

    def initPreweight(self, pathPreWeight: str = None, rmModule=True):

        preW = self.getWeight(pathPreWeight)
        assert preW is not None, 'weighth in {} is empty'.format(pathPreWeight)
        modelW = self.state_dict()
        preWDict = OrderedDict()
        # modelWDict = OrderedDict()

        for k, v in preW.items():
            if rmModule:
                preWDict[k.replace('module.', "")] = v
            else:
                preWDict[k] = v

        shareW = {k: v for k, v in preWDict.items() if str(k) in modelW}
        assert shareW, 'shareW is empty'
        self.load_state_dict(preWDict, strict=False)

    @staticmethod
    def getWeight(pathPreWeight: str = None):
        if pathPreWeight is not None:
            return torch.load(pathPreWeight, map_location=torch.device('cpu'))
        else:
            return None

    @staticmethod
    def padToScale(img, netScale):
        _, _, h, w = img.size()
        oh = int(np.ceil(h * 1.0 / netScale) * netScale)
        ow = int(np.ceil(w * 1.0 / netScale) * netScale)
        img = F.pad(img, [0, ow - w, 0, oh - h], mode='reflect')
        return img

class conv3x3(nn.Module):
    def __init__(self, in_c, out_c, cfg, ks=2):
        super(conv3x3, self).__init__()

        self.conv = nn.Sequential(
            nn.ReplicationPad2d([1, 1, 1, 1]),
            nn.Conv2d(in_channels=in_c, out_channels=out_c, kernel_size=3, stride=ks, bias=False)
        )
        self.norm2d = getNorm(cfg.netNorm, num_features=out_c)
        self.activate = getAct(cfg.netActivate, num_parameters=out_c)

    def forward(self, x):
        # x = F.pad(input=x, pad=[2, 2, 2, 2], mode='replicate')
        x = self.conv(x)
        x = self.norm2d(x)
        x = self.activate(x)

        return x

class deconv(nn.Module):
    def __init__(self, in_c, out_c, kernel_size, cfg=None):
        super(deconv, self).__init__()

        if isinstance(kernel_size, tuple):
            kernel_size = kernel_size
        else:
            kernel_size = (kernel_size, kernel_size)
        self.bias = False if 'identity' == cfg.netNorm else True

        self.padList = [(kernel_size[1] - 1) // 2, (kernel_size[1] - 1) // 2,
                        (kernel_size[0] - 1) // 2, (kernel_size[0] - 1) // 2]

        self.deconv = nn.Sequential(
            nn.UpsamplingBilinear2d(scale_factor=2),
            nn.ReplicationPad2d(self.padList),
            nn.Conv2d(in_channels=in_c, out_channels=out_c, kernel_size=kernel_size, stride=1, bias=False)
        )

        self.norm2d = getNorm(cfg.netNorm, num_features=out_c)
        self.activate = getAct(cfg.netActivate, num_parameters=out_c)

    def forward(self, x: torch.Tensor, skip: torch.Tensor):

        x = self.deconv(x)
        x = self.norm2d(x)
        x = self.activate(x)

        return torch.cat([x, skip], dim=1)

class FuseLayer(nn.Module):
    def __init__(self, c_h, c_ef):
        super(FuseLayer, self).__init__()
        self.c_ef = c_ef
        self.c_h = c_h

        self.norm = nn.InstanceNorm2d(c_h, affine=False)

        self.convE1 = nn.Conv2d(c_ef, c_h, kernel_size=1, stride=1, padding=0, bias=True)
        self.convE2 = nn.Conv2d(c_ef, c_h, kernel_size=1, stride=1, padding=0, bias=True)

        self.convF1 = nn.Conv2d(c_ef, c_h, kernel_size=1, stride=1, padding=0, bias=True)
        self.convF2 = nn.Conv2d(c_ef, c_h, kernel_size=1, stride=1, padding=0, bias=True)

        self.convMask = nn.Sequential(nn.ReplicationPad2d(padding=[1, 1, 1, 1]),
                                      nn.Conv2d(c_h, 1, kernel_size=3, stride=1, bias=True),
                                      nn.Sigmoid())

    def forward(self, h_in, z_e, z_f):
        h = self.norm(h_in)

        gammaE = self.convE1(z_e)
        betaE = self.convE2(z_e)
        E = gammaE * h + betaE

        gammaF = self.convF1(z_f)
        betaF = self.convF2(z_f)
        F = gammaF * h + betaF

        hEF = h
        M = self.convMask(hEF)
        out = M * E + (1.0 - M) * F

        return out

class FuseBlock(nn.Module):
    def __init__(self, cin, cout, c_ef):
        super(FuseBlock, self).__init__()
        self.cin = cin
        self.cout = cout

        self.AAD1 = FuseLayer(cin, c_ef)
        self.conv1 = nn.Sequential(
            nn.LeakyReLU(negative_slope=0.2, inplace=True),
            nn.ReplicationPad2d(padding=[1, 1, 1, 1]),
            nn.Conv2d(cin, cin, kernel_size=3, stride=1, bias=False)
        )

        self.AAD2 = FuseLayer(cin, c_ef)
        self.conv2 = nn.Sequential(
            nn.LeakyReLU(negative_slope=0.2, inplace=True),
            nn.ReplicationPad2d(padding=[1, 1, 1, 1]),
            nn.Conv2d(cin, cout, kernel_size=3, stride=1, bias=False)
        )

        if cin != cout:
            self.AAD3 = FuseLayer(cin, c_ef)
            self.conv3 = nn.Sequential(
                nn.LeakyReLU(negative_slope=0.2, inplace=True),
                nn.ReplicationPad2d(padding=[1, 1, 1, 1]),
                nn.Conv2d(cin, cout, kernel_size=3, stride=1, bias=False)
            )

    def forward(self, h, z_e, z_f):

        x = self.AAD1(h, z_e, z_f)
        x = self.conv1(x)

        x = self.AAD2(x, z_e, z_f)
        x = self.conv2(x)

        if self.cin != self.cout:
            h = self.AAD3(h, z_e, z_f)
            h = self.conv3(h)
        x = x + h

        return x

class PreActBlock(nn.Module):
    # Pre-activation version of the BasicBlock

    def __init__(self, inC, outC, cfg, stride=1):
        super(PreActBlock, self).__init__()

        # self.norm2d = getNorm(cfg.netNorm, num_features=out_c)
        # self.activate = getAct(cfg.netActivate, num_parameters=out_c)
        self.expansion = 1
        self.norm1 = nn.GroupNorm(4, inC, affine=True)
        self.conv1 = nn.Conv2d(inC, outC, kernel_size=3, stride=stride, padding=1, bias=False)
        self.norm2 = nn.GroupNorm(4, outC, affine=True)
        self.conv2 = nn.Conv2d(outC, outC, kernel_size=3, stride=1, padding=1, bias=False)
        self.activate = nn.LeakyReLU(0.2, True)
        if stride != 1 or inC != self.expansion * outC:
            self.shortcut = nn.Sequential(
                nn.Conv2d(inC, self.expansion * outC, kernel_size=1, stride=stride, bias=False)
            )

    def forward(self, x):
        out = self.activate(self.norm1(x))
        shortcut = self.shortcut(out) if hasattr(self, 'shortcut') else x
        out = self.conv1(out)
        out = self.conv2(self.activate(self.norm2(out)))
        out += shortcut
        return out

class PreActBottleneck(nn.Module):
    # Pre-activation version of the original Bottleneck module.

    def __init__(self, inC, outC, stride=1):
        super(PreActBottleneck, self).__init__()
        midC = inC // 2 if outC == inC else inC
        self.norm1 = nn.GroupNorm(4, inC, affine=True)
        self.conv1 = nn.Conv2d(inC, midC, kernel_size=1, bias=False)

        self.norm2 = nn.GroupNorm(4, midC, affine=True)
        self.conv2 = nn.Conv2d(midC, midC, kernel_size=3, stride=stride, padding=1, bias=False)

        self.norm3 = nn.GroupNorm(4, midC, affine=True)
        self.conv3 = nn.Conv2d(midC, outC, kernel_size=1, bias=False)

        self.activate = nn.LeakyReLU(0.2, True)

        if stride != 1 or inC != outC:
            self.shortcut = nn.Sequential(
                nn.Conv2d(inC, outC, kernel_size=1, stride=stride, bias=False)
            )

    def forward(self, x):
        out = self.activate(self.norm1(x))
        shortcut = self.shortcut(out) if hasattr(self, 'shortcut') else x
        out = self.conv1(out)
        out = self.conv2(self.activate(self.norm2(out)))
        out = self.conv3(self.activate(self.norm3(out)))
        out += shortcut
        return out

# end fuseNet --------------------------------------------------------------------------------------

# for PostProcess ---------------------------------------------------------------------------------

def conv1x1(in_channels, out_channels, stride=1):
    return nn.Conv2d(in_channels, out_channels, kernel_size=1, stride=stride, padding=0, bias=True)

def conv2D3x3(in_channels, out_channels, stride=1):
    return nn.Conv2d(in_channels, out_channels, kernel_size=3, stride=stride, padding=1, bias=True)

class ResBlock(nn.Module):
    def __init__(self, in_channels, out_channels, stride=1, downsample=None, res_scale=1):
        super(ResBlock, self).__init__()
        self.res_scale = res_scale
        self.conv1 = conv2D3x3(in_channels, out_channels, stride)
        self.relu = nn.ReLU(inplace=True)
        self.conv2 = conv2D3x3(out_channels, out_channels)

    def forward(self, x):
        x1 = x
        out = self.conv1(x)
        out = self.relu(out)
        out = self.conv2(out)
        out = out * self.res_scale + x1
        return out

# end PostProcess ---------------------------------------------------------------------------------

if __name__ == '__main__':
    from configs import configEVI

    cfg = configEVI.Config()
    pass

max_stars_repo_name: YU-Zhiyang/WEVI

max_issues_repo_name: YU-Zhiyang/WEVI

max_forks_repo_name: YU-Zhiyang/WEVI


text:
# -*- coding: utf-8 -*-
#find the integral of a function f(x) from a to b i.e.

# ∫from a to b f(x)dx

#In python we use numerical quadrature to achieve this with the
# scipy.integrate.quad command. As a specific example, lets integrate

# y=x(squared)

#from x=0 to x=1.
#You should be able to work out that the answer is 1/3.

from scipy.integrate import quad

def integrand(x):
    return x**2

ans, err = quad(integrand, 0, 1)
print ans

max_stars_repo_name: lastralab/Statistics

max_issues_repo_name: lastralab/Statistics

max_forks_repo_name: lastralab/Statistics

reason: from scipy


text:
from summit.strategies import TSEMO
from summit.benchmarks import DTLZ2
from summit.utils.models import GPyModel
from GPy.kern import Exponential
import numpy as np
import json
import warnings
from fastprogress.fastprogress import master_bar, progress_bar

#Two tests were run; one with 1500 and one with 4000
N_SPECTRAL_POINTS = 4000
DATE="20200618"
tsemo_options = dict(pop_size=100,                          #population size for NSGAII
                     iterations=100,                        #iterations for NSGAII
                     n_spectral_points=N_SPECTRAL_POINTS,   #number of spectral points for spectral sampling
                     num_restarts=200,                      #number of restarts for GP optimizer (LBSG)
                     parallel=True)                         #operate GP optimizer in parallel
description="Description: Used exponential kernel instead of matern and increase num restarts. Also, catch and save errors."

def dtlz2_test():
    #Run the DTLZ2 benchmark
    errors = 0
    num_inputs=6
    num_objectives=2
    lab = DTLZ2(num_inputs=num_inputs,
                num_objectives=num_objectives)
    models = {f'y_{i}': GPyModel(Exponential(input_dim=num_inputs,ARD=True))
              for i in range(num_objectives)}

    warnings.filterwarnings("ignore", category=RuntimeWarning)
    tsemo = TSEMO(lab.domain, models=models, random_rate=0.00)
    experiments = tsemo.suggest_experiments(5*num_inputs)

    mb = master_bar(range(1))
    for j in mb:
        mb.main_bar.comment = f'Repeats'
        for i in progress_bar(range(100), parent=mb):
            mb.child.comment = f'Iteration'
            # Run experiments
            experiments = lab.run_experiments(experiments)

            # Get suggestions
            try:
                experiments = tsemo.suggest_experiments(1, experiments,
                                                        **tsemo_options)
            except Exception as e:
                print(e)
                errors +=1

        tsemo.save(f'new_tsemo_params_{j}.json')

if __name__ == '__main__':
    tsemo_options.update({'description': description})
    with open(f'new_params.json', 'w') as f:
        json.dump(tsemo_options,f)
    dtlz2_test()

max_stars_repo_name: dswigh/summit

max_issues_repo_name: dswigh/summit

max_forks_repo_name: dswigh/summit


text:
import numpy as np

from sklearn.base import BaseEstimator, TransformerMixin
from typing import Sequence, Optional, Union, Tuple

class Segmented:
    def __init__(self,
                 breaks: Union[int, Sequence] = (),
                 break_method: Optional[str] = None):
        self.breaks = breaks
        self.break_method = break_method
        self.break_locations_ = None

    def _choose_breaks(self,
                       X: np.ndarray,
                       y: Optional[np.ndarray],
                       num: int) -> Tuple[float, ...]:
        break_method = (self.break_method or 'quantile').lower()

        if break_method == 'quantile':
            q = np.linspace(0, 1, num + 2)[1:-1]
            breaks = set(np.nanquantile(X, q))
            if len(breaks) != num:
                raise RuntimeError("quantile method returned duplicate breaks")
            return tuple(sorted(breaks))
        else:
            raise NotImplementedError

    def fit(self, X: np.ndarray, y: Optional[np.ndarray] = None):
        X = self._validate_X(X)
        assert len(np.unique(X)) > 1

        if isinstance(self.breaks, (list, tuple, np.ndarray)):
            assert len(self.breaks) == len(set(self.breaks))
            self.break_locations_ = tuple(self.breaks)
        else:
            self.break_locations_ = self._choose_breaks(X, y, self.breaks)

        return self

    def _validate_X(self, X: np.ndarray) -> np.ndarray:
        if hasattr(X, 'values'):
            X = X.values
        if len(X.shape) == 2:
            if X.shape[1] != 1:
                raise RuntimeError("Expected X to be vector")
            X = X.squeeze(-1)
        return X

class SegmentedBoundary(Segmented):
    def __init__(self,
                 breaks: Union[int, Sequence] = (),
                 boundary_breaks: Optional[Sequence] = None,
                 break_method: Optional[str] = None):
        self.boundary_breaks = boundary_breaks
        super().__init__(breaks=breaks, break_method=break_method)

    def fit(self, X: np.ndarray, y: Optional[np.ndarray] = None):
        super().fit(X, y)
        s, e = self._choose_boundary(X)

        if not all(s < k < e for k in self.break_locations_):
            raise ValueError(f"Some breaks fall outside boundary: {s, e}")

        self.break_locations_ = (s,) + self.break_locations_ + (e,)

        return self

    def _choose_boundary(self, X: np.ndarray) -> Tuple:
        if isinstance(self.boundary_breaks, (list, tuple, np.ndarray)):
            return self.boundary_breaks
        else:
            return np.nanmin(X), np.nanmax(X)

class PiecewiseCumulative(Segmented, BaseEstimator, TransformerMixin):
    def transform(self, X: np.ndarray, y=None):
        X = self._validate_X(X)
        out = [X]
        raise RuntimeError("TODO: why [1:]?")
        for brk in self.break_locations_[1:]:
            out.append(np.maximum(X - brk, 0))

        out = np.vstack(out).T
        return out

class Piecewise(SegmentedBoundary, BaseEstimator, TransformerMixin):
    def __init__(self,
                 breaks: Union[int, Sequence] = (),
                 boundary_breaks: Optional[Sequence] = None,
                 break_method: Optional[str] = None,
                 dampen: float = 0.0):
        """

        :param breaks: Either (a) integer for number of breaks to find, if break_method supports this, (b) actual
        breakpoints. The piecewise function has num-breaks + 1 degrees of freedom.
        :param break_method: Method for finding breaks; default is to use quantiles. Ignored if actual breakpoints are
        specified w/`breaks` argument.
        :param boundary_breaks: In addition to the breaks chosen with `breaks`, these outer breaks will be added. The
        function's behavior outside of these breaks is determined by the `dampen` argument. TODO: If dampen is 0.0,
        placement of boundary breaks has no effect? TODO: default behavior is minmax is dampen is zero, else (TODO)
        :param dampen: A float between 0 and 1 inclusive. If 1.0, then the function is flat outside of the boundary-
        breaks. If 0.0, then the function extrapolates the slope from before the boundary.
        :param control_break_idx: TODO
        """
        assert 0.0 <= dampen <= 1.0
        self.dampen = dampen
        super().__init__(breaks=breaks, break_method=break_method, boundary_breaks=boundary_breaks)

    def fit(self, X: np.ndarray, y: Optional[np.ndarray] = None):
        super().fit(X, y)

    def transform(self, X: np.ndarray, y=None):
        X = self._validate_X(X)

        num_knots = len(self.break_locations_)
        if num_knots == 2:
            return X

        spans = np.diff(self.break_locations_)
        spans /= spans.mean()

        out = []

        # LHS:
        for i in range(1, middle_idx + 1):
            clipped = np.clip(X, None, self.break_locations_[i])
            x = self._dampen(clipped, 'lhs') - self.break_locations_[i]
            if self.scale_by_span:
                x *= spans[i - 1]
            out.append(x)

        # RHS:
        for i in range(middle_idx, num_knots - 1):
            clipped = np.clip(X, self.break_locations_[i], None)
            x = self._dampen(clipped, 'rhs') - self.break_locations_[i]
            if self.scale_by_span:
                x *= spans[i]
            out.append(x)

        return np.vstack(out).T

    def _dampen(self, clipped: np.ndarray, side: str) -> np.ndarray:
        # TODO: could make dampening increase with distance?
        x = clipped.copy()
        w1 = self.dampen
        w2 = 1 - w1
        if side == 'rhs':
            outer_break = self.break_locations_[-1]
            is_outer = x > outer_break
        else:
            outer_break = self.break_locations_[0]
            is_outer = x < outer_break
        x[is_outer] = (w1 * outer_break + w2 * x)[is_outer]
        return x

    def _choose_boundary(self, X: np.ndarray) -> Tuple:
        s, e = super()._choose_boundary(X)
        if not len(self.break_locations_) or self.dampen == 0.0:
            return s, e
        raise RuntimeError("TODO: mean of (last_quantile,100.0), etc.")
        return s, e

try:
    from patsy.mgcv_cubic_splines import _get_free_crs_dmatrix, _get_centering_constraint_from_dmatrix, _get_crs_dmatrix
except ImportError:
    def _no_patsy(*args, **kwargs):
        raise ImportError("Must install `patsy` for cardinal basis")

    _get_free_crs_dmatrix = _get_centering_constraint_from_dmatrix = _get_crs_dmatrix = _no_patsy

class NaturalSpline(SegmentedBoundary, BaseEstimator, TransformerMixin):
    def __init__(self,
                 knots: Union[int, Sequence] = (),
                 boundary_knots: Optional[Sequence] = None,
                 knot_select_method: Optional[str] = None,
                 basis: str = 'cardinal'):
        self.basis = basis
        self.constraints_ = None
        super().__init__(breaks=knots, boundary_breaks=boundary_knots, break_method=knot_select_method)

    def fit(self, X: np.ndarray, y: Optional[np.ndarray] = None):
        super().fit(X, y)
        if self.basis.lower() == 'cardinal':
            self.constraints_ = _get_centering_constraint_from_dmatrix(_get_free_crs_dmatrix(x, self.break_locations_))
        return self

    def transform(self, X: np.ndarray, y=None):
        X = self._validate_X(X)
        if self.basis.lower().startswith('trunc'):
            return self._trunc_power_basis(X)
        elif self.basis.lower() == 'cardinal':
            return self._cardinal_basis(X)
        else:
            raise ValueError(f"Unrecognized basis: {self.basis}")

    def _cardinal_basis(self, X) -> np.ndarray:
        return _get_crs_dmatrix(X, knots=self.break_locations_, constraints=self.constraints_)

    def _trunc_power_basis(self, X: np.ndarray) -> Sequence:
        out = [X]
        sub = self._d(X, self.break_locations_[-2])
        for knot in self.break_locations_[:-2]:
            out.append(self._d(X, knot) - sub)
        return np.vstack(out).T

    def _d(self, X: np.ndarray, knot: float) -> np.ndarray:
        last_knot = self.break_locations_[-1]
        p1 = np.maximum(X - knot, 0) ** 3
        p2 = np.maximum(X - last_knot, 0) ** 3
        return (p1 - p2) / (last_knot - knot)

max_stars_repo_name: strongio/strong-glm

max_issues_repo_name: strongio/strong-glm

max_forks_repo_name: strongio/strong-glm