text:
#!/usr/bin/env python
# -*- coding: utf-8 -*-

"""
Creating violin plots in PGFplots (two-sided version)

Based on:
https://matplotlib.org/_modules/matplotlib/axes/_axes.html#Axes.violinplot
https://github.com/statsmodels/statsmodels/blob/master/statsmodels/graphics/boxplots.py

Author: Gertjan van den Burg

"""

import argparse
import json
import numpy as np
import os
import math

from scipy.stats import gaussian_kde

from .core import (
    CORPUS_NAMES,
    ORDERED_DETECTORS,
    check_detectors,
    clean_detector_name,
)
from .latex import build_latex_doc

# Color
# COLOR_LEFT = "B40204"
# COLOR_RIGHT = "00AABA"
# COLOR_MINMAX = "FF0000"

# Grayscale
COLOR_LEFT = "101010"
COLOR_RIGHT = "878787"
COLOR_MINMAX = "000000"

USE_LOG = True

def transform(x):
    if USE_LOG:
        return math.log(x, 10)
    return x

def untransform(x):
    if USE_LOG:
        return pow(10, x)
    return x

def _interpolate(coords, values, x):
    """ Return the estimated value of x by interpolating the nearest neighbors
    in coords. It is assumed coords is sorted and is of the same length as
    values.
    """
    if x in coords:
        return values[coords == x]
    below_idx, above_idx = None, None
    for idx, c in enumerate(coords):
        if c < x:
            below_idx = idx
        if c > x:
            above_idx = idx
            break
    avg_val = (values[below_idx] + values[above_idx]) / 2
    return avg_val

def _single_violin_data(pos, pos_data, width, side, plot_opts):
    # Based almost entirely on_single_violin from statsmodels
    bw_factor = plot_opts.get("bw_factor", None)

    def _violin_range(pos_data, plot_opts):
        """Return array with correct range, with which violins can be plotted."""
        cutoff = plot_opts.get("cutoff", False)
        cutoff_type = plot_opts.get("cutoff_type", "std")
        cutoff_val = plot_opts.get("cutoff_val", 1.5)

        s = 0.0
        if not cutoff:
            if cutoff_type == "std":
                s = cutoff_val * np.std(pos_data)
            else:
                s = cutoff_val

        x_lower = kde.dataset.min() - s
        x_upper = kde.dataset.max() + s
        return np.linspace(x_lower, x_upper, 501)

    pos_data = np.asarray(pos_data)
    kde = gaussian_kde(pos_data, bw_method=bw_factor)

    xvals = _violin_range(pos_data, plot_opts)
    violin = kde.evaluate(xvals)

    # NOTE: we removed normalization by violin.max()
    violin = width * violin

    if side == "both":
        envelope_l, envelope_r = (-violin + pos, violin + pos)
    elif side == "right":
        envelope_l, envelope_r = (np.zeros_like(violin) + pos, violin + pos)
    elif side == "left":
        envelope_l, envelope_r = (-violin + pos, np.zeros_like(violin) + pos)
    else:
        msg = "`side` parameter should be one of {'left', 'right', 'both'}."
        raise ValueError(msg)

    return xvals, envelope_l, envelope_r

def get_median_coords(coords, left, right, median):
    data = {}
    data["xleft"] = _interpolate(coords, left, median)
    data["xright"] = _interpolate(coords, right, median)
    data["yleft"] = median
    data["yright"] = median
    return data

def get_extrema_coords(pos, pos_data, width, side):
    # min
    xleft = pos
    xleft -= width if side in ["left", "both"] else 0
    xright = pos
    xright += width if side in ["right", "both"] else 0
    yleft = yright = np.min(pos_data)
    min_coords = {
        "xleft": xleft,
        "xright": xright,
        "yleft": yleft,
        "yright": yright,
    }
    # max
    yleft = yright = np.max(pos_data)
    max_coords = {
        "xleft": xleft,
        "xright": xright,
        "yleft": yleft,
        "yright": yright,
    }
    return min_coords, max_coords

def generate_violin_data(
    summary_data, side="both", showmedian=True, showextrema=True, plot_opts={}
):

    check_detectors(summary_data["runtimes"].keys())

    dataset = list(
        map(
            np.asarray,
            [
                list(map(transform, summary_data["runtimes"][key]))
                for key in ORDERED_DETECTORS
            ],
        )
    )

    positions = np.arange(len(dataset)) + 1
    pos_span = np.max(positions) - np.min(positions)
    width = np.min(
        [0.15 * np.max([pos_span, 1.]), plot_opts.get("violin_width", 0.8) / 2.]
    )

    violin_data = []
    for pos_data, pos, name in zip(dataset, positions, ORDERED_DETECTORS):
        xvals, envelope_l, envelope_r = _single_violin_data(
            pos, pos_data, width, side, plot_opts
        )

        # return back to actual data
        xvals = np.array([untransform(x) for x in xvals])
        pos_data = np.array([untransform(x) for x in pos_data])

        data = {
            "name": name,
            "side": side,
            "xvals": xvals,
            "envelope_l": envelope_l,
            "envelope_r": envelope_r,
        }

        if showmedian:
            data["median"] = get_median_coords(
                xvals, envelope_l, envelope_r, np.median(pos_data)
            )
        if showextrema:
            data["min"], data["max"] = get_extrema_coords(
                pos, pos_data, width / 3, side
            )

        violin_data.append(data)

    return violin_data

def generate_tex_for_line(xleft=0, yleft=0, xright=0, yright=0, linestyle=""):
    tex = ""
    tex += "\\addplot[%s] coordinates {%%\n" % linestyle
    tex += "(%.16f, %.16f)\n" % (xleft, yleft)
    tex += "(%.16f, %.16f)\n" % (xright, yright)
    tex += "};\n"
    return tex

def generate_tex_for_violin(
    violin, edgecolor=None, edgethick=None, fillcolor=None, alpha=0.5
):
    name = violin["name"] + violin["side"]

    edgecolor = "none" if edgecolor is None else edgecolor
    edgethick = "" if edgethick is None else ", " + edgethick
    fillcolor = "fill=none" if fillcolor is None else fillcolor
    left_name, right_name = name + "Left", name + "Right"

    tex = "\\addplot [draw=%s %s, name path=%s] coordinates {%%\n" % (
        edgecolor,
        edgethick,
        left_name,
    )
    for xx, yy in zip(violin["envelope_l"], violin["xvals"]):
        tex += "(%.16f, %.16f)\n" % (xx, yy)
    tex += "};\n"
    tex += "\\addplot [draw=%s %s, name path=%s] coordinates {%%\n" % (
        edgecolor,
        edgethick,
        right_name,
    )
    for xx, yy in zip(violin["envelope_r"], violin["xvals"]):
        tex += "(%.16f, %.16f)\n" % (xx, yy)
    tex += "};\n"
    tex += "\\addplot [%s, opacity=%f] fill between [of=%s and %s];\n" % (
        fillcolor,
        alpha,
        left_name,
        right_name,
    )

    if "median" in violin:
        # linestyle = "dashed, dash pattern=on 2pt off 2pt"
        violin["median"]["linestyle"] = "densely dotted, thick, black"

        tex += generate_tex_for_line(**violin["median"])
    if "min" in violin:
        violin["min"]["linestyle"] = "solid, ColorMinMax"
        tex += generate_tex_for_line(**violin["min"])
    if "max" in violin:
        violin["max"]["linestyle"] = "solid, ColorMinMax"
        tex += generate_tex_for_line(**violin["max"])

    return tex

def generate_latex(violindata, legend_data, opacity=0.5):
    abbrev = [clean_detector_name(d) for d in ORDERED_DETECTORS]
    xtick = ",".join([str(i + 1) for i in range(len(abbrev))])
    xticklabels = ",".join(abbrev)

    yrange = [pow(10, x) for x in [-6, -4, -2, 0, 2, 4]]
    ytick = ",".join([str(i) for i in yrange])

    legend_entries = ", ".join(
        [CORPUS_NAMES.get(c) for c in legend_data["corpora"]]
    )

    tex = (
        "\\documentclass[preview=true]{standalone}\n"
        "\\pdfinfoomitdate=1\n"
        "\\pdftrailerid{}\n"
        "\\pdfsuppressptexinfo=1\n"
        "\\usepackage{tikz}\n"
        "\\usepackage{pgfplots}\n"
        "\\pgfplotsset{compat=1.16}\n"
        "\\usepgfplotslibrary{fillbetween}\n"
        "\\definecolor{ColorLeft}{HTML}{%s}\n"
        "\\definecolor{ColorRight}{HTML}{%s}\n"
        "\\definecolor{ColorMinMax}{HTML}{%s}\n"
        "\\begin{document}\n"
        "\\begin{tikzpicture}\n"
        "\\begin{semilogyaxis}[\n"
        "xtick={%s},\n"
        "xticklabels={%s},\n"
        "ytick={%s},\n"
        "ylabel={Runtime (s)},\n"
        "width=600pt,\n"
        "height=200pt,\n"
        "ymajorgrids,\n"
        "grid style={opacity=0.1},\n"
        "legend entries={%s},\n"
        "legend pos={south west},\n"
        "]\n"
        % (
            COLOR_LEFT,
            COLOR_RIGHT,
            COLOR_MINMAX,
            xtick,
            xticklabels,
            ytick,
            legend_entries,
        )
    )

    tex += (
        "\\addlegendimage{only marks, mark=square*, ColorLeft, opacity=%g}\n"
        % (opacity)
    )
    tex += (
        "\\addlegendimage{only marks, mark=square*, ColorRight, opacity=%g}\n"
        % (opacity)
    )

    for corpus in violindata:
        for violin in violindata[corpus]:
            fillcolor = (
                "ColorLeft" if violin["side"] == "left" else "ColorRight"
            )
            tex += generate_tex_for_violin(
                violin, edgecolor="black", fillcolor=fillcolor, alpha=opacity
            )

    tex += "\\end{semilogyaxis}\n" "\\end{tikzpicture}\n" "\\end{document}"
    return tex

def create_twosided_violin(corpus_data, output_file):
    corpora = sorted(corpus_data.keys())
    sides = ["left", "right"]
    assert len(corpora) == 2
    legend_data = {"corpora": corpora, "colors": [COLOR_LEFT, COLOR_RIGHT]}

    violindata = {}
    for corpus, side in zip(corpora, sides):
        violindata[corpus] = generate_violin_data(
            corpus_data[corpus], side=side, showmedian=True
        )

    tex = generate_latex(violindata, legend_data)
    tex_file = os.path.splitext(output_file)[0] + ".tex"
    with open(tex_file, "w") as fid:
        fid.write(tex)
    build_latex_doc(tex, output_name=output_file)

def parse_args():
    parser = argparse.ArgumentParser()
    parser.add_argument(
        "-o", dest="output", help="Output pdf file to write to", required=True
    )
    parser.add_argument(
        "-s",
        dest="summaries",
        help="Summary file(s) with the results",
        required=True,
        nargs="+",
    )
    return parser.parse_args()

def main():
    args = parse_args()
    all_data = {}
    for summary_file in args.summaries:
        with open(summary_file, "r") as fid:
            data = json.load(fid)
        all_data[data["corpus"]] = data

    create_twosided_violin(all_data, args.output)


text:
"""
    Extract features from DBSCAN clusters that show results related to the quality of clusters and number of images
	author: Cristina Luna, Ricardo Kleinlein
	date: 05/2020
	Usage:
		python3 ClusterScan.py
		PARAMETERS FOR RUNNING
         --root-input-dir ../supervised_results/BASELINE_OFICIAL_LN24H_SOLOIND
         --add-noise-class False
         --assignation-type x
         --qualityTh-valid-cluster 0.0
         --assignation-method numImgs
         --dataset Google
         --program-participants-folder ../../data/datasets/DATASET_GOOGLE_IMGS/participants
         --output-dir ../../data/datasets/popularity_module_features/aux
         --with-previous-individual-clusters
         --first-param-name eps
         --second-param-name min_samples
         --quality-metric silhouette
         --root-path-MTCNN-results ../DATASET_GOOGLE_IMGS/VIDEO_DB_MTCNN
         --program-name TOTAL
	Options:
	    --root-path-MTCNN-results: Root path to the results obained by MTCNN (embs, bbox, mtcnn_debub, bbox_summ)
	    --program-participants-folder: Path to the folder with the .csv that contain the names of the
	    participants per program (See /AQuery_launcher/Aqeury_launcher_bingEngine.py
	    --first-param-name: Name of the first param as it appears in the DBSCAN and HDBSCAN libraries.
                 -DBSCAN: https://scikit-learn.org/stable/modules/generated/sklearn.cluster.DBSCAN.html
                 -HDBSCAN: https://github.com/scikit-learn-contrib/hdbscan
                 [default: eps (DBSCAN param.)]
	    --first-param-lower: Lower value for first parameter [default: 0.5]. For DBSCAN this parameter is eps
	    --first-param-upper: Upper value for fisrt parameter [default: 12]. For DBSCAN this parameter is eps
	    --first-param-step: Step between 2 consecutive values of the first parameter. For DBSCAN this parameter is eps
	    --second-param-value: List with second param values [default: [10]]
	    --metric: Distance metric in DBSCAN [Options: euclidean, cosine, manhattan, l1, l2]
	    --dataset: Dataset reference name to extract clustering & save results(e.g Google, OCR, Google_OCR, Google_Program ...)
		--output-dir	Directory to save results in.By default, it will act as root folder and results will be saved in the
		path: output_dir/Cfacial_clustering_results/'dataset'/
		                                                |_cluster_models/'program'/
		                                                                    |_cluster_eps_x_minSamples_x.sav
		                                                |_cluster_graphics/parameters/'program'
		                                                                              |_(matrix/graphics...)
		--individual-clusters: True if we want a cluster per participant, else it creates a common cluster for all participants
        --quality-metric: Metric to use as quality and select the best cluster combination of parameters. Options ( cluster, v_meassure, silhouette, calinski_harabasz).
        --with-previous-individual-clusters: True if we did a previous individual cluster per participant, else False
        --with-previous-common-clusters: True if we did a previous common cluster with all the participants, else False.
		-h, --help	Display script additional help
"""

import os
import pandas as pd
import numpy as np
import src.metrics.metrics_clustering as metrics
import sklearn
from src.BaseArgs import ClusteringCleanerArgs
import src.utils.loader as loader

class ClusterFeatures():
    def __init__(self, input_dir, root_path_MTCNN_results, output_dir,program_name,program_participants, root_input_dir_parameters,
                 input_path_models,input_path_parameters, first_param_name = "eps", second_param_name ="min_samples", quality_metric="silhouette"):
        self.program = program_name
        self.program_participants = program_participants
        self.root_input_dir_parameters = root_input_dir_parameters
        self.quality_metric = quality_metric
        self.input_dir = input_dir
        self.first_param_name = first_param_name
        self.second_param_name = second_param_name
        self.root_input_path_embs = os.path.join(root_path_MTCNN_results, program_name, "embeddings_sum")
        self.root_input_path_bbox = os.path.join(root_path_MTCNN_results, program_name, "boundingboxes_sum")
        self.root_input_path_mtcnn = os.path.join(root_path_MTCNN_results, program_name, "mtcnn_debug_sum")
        self.input_path_models = input_path_models
        self.input_path_parameters = input_path_parameters

    def extract_clt_features(self,participant):
        """
        Extract DBSCAN cluster features
        :param participant: Name of the participant in folder
        :return: 1. dataframe with statistics of global cluster features and
                 2. dataframe with statistics of individual cluster features
        """
        #Load best parameters per participant
        input_path_parameters_participant = os.path.join(self.input_path_parameters, participant)
        input_path_models_participant = os.path.join(self.input_path_models, participant)
        _, _, clt_embs, embs_total_labels = loader.load_bboxes_and_embs(
            self.root_input_path_bbox, self.root_input_path_embs, [participant])
        df_global_features = pd.DataFrame([], columns=["id", "num_clts", "avg_silhouette", "avg_intra_clt_dist", "avg_inter_clt_dist", "faces_in_noise_clt", "num_core_samples"])
        df_individual_clt_features = pd.DataFrame([], columns=["id", "clt_id", "n_faces",
                               "avg_silhouette","std_silhouette",
                               "avg_intra_clt_dist", "std_intra_clt_dist",
                               "avg_inter_clt_dist", "std_inter_clt_dist",
                               "num_core_samples"])
        try:
            first_param,snd_param = loader.load_best_parameters_new(input_path_parameters_participant, self.quality_metric)
            #Load cluster
            path_cluster = os.path.join(input_path_models_participant,"cluster_"+self.first_param_name+"_"+str(first_param)+"_"+
                                        self.second_param_name+"_"+str(snd_param)+".sav")
            clt = loader.load_cluster(path_cluster)
            clt_labels = clt.labels_
            #clt_core_samples = clt.components_
            n_core_samples = clt.core_sample_indices_
            labels_core_samples = clt_labels[clt.core_sample_indices_]
            labels_ID = list(set(np.unique(clt_labels)).difference([-1]))
            intra_clust_dists, inter_clust_dists = metrics.intra_cluster_distance(clt_embs, clt_labels)
            silh_samples = sklearn.metrics.silhouette_samples(clt_embs, clt_labels)
            #GLOBAL FEATURES:
            faces_in_noise_clt = np.where(clt_labels==-1)
            silh_global = np.nanmean(silh_samples)
            avg_intra_clt_dist = np.nanmean(intra_clust_dists)
            avg_inter_clt_dist = np.nanmean(inter_clust_dists)
            df_global_features = pd.DataFrame([[participant, len(labels_ID),silh_global,avg_intra_clt_dist,avg_inter_clt_dist,len(faces_in_noise_clt[0]),len(n_core_samples)]], columns=["id", "num_clts", "avg_silhouette", "avg_intra_clt_dist", "avg_inter_clt_dist", "faces_in_noise_clt", "num_core_samples"])

            # PER CLUSTER FEATURES:
            for label in labels_ID:
                index_current_label = np.where(clt_labels==label)
                df_individual_clt_features = df_individual_clt_features.append(pd.DataFrame([[participant, label, len(index_current_label[0]),
                   np.nanmean(silh_samples[index_current_label]),np.nanstd(silh_samples[index_current_label]),
                   np.nanmean(intra_clust_dists[index_current_label]),np.nanstd(intra_clust_dists[index_current_label]),
                   np.nanmean(inter_clust_dists[index_current_label]),np.nanstd(inter_clust_dists[index_current_label]),
                   len(np.where(labels_core_samples == label)[0])]],
                   columns=["id", "clt_id", "n_faces",
                               "avg_silhouette","std_silhouette",
                               "avg_intra_clt_dist", "std_intra_clt_dist",
                               "avg_inter_clt_dist", "std_inter_clt_dist",
                               "num_core_samples"]))
        except FileNotFoundError:
            print("PASS")
            pass
        return df_global_features, df_individual_clt_features

if __name__ == "__main__":
    #Parse arguments
    clustering_cleaner_args_obj = ClusteringCleanerArgs()
    args = clustering_cleaner_args_obj.parse()

    # Get input path
    if (args.assignation_type == "x"):
        args.qualityTh_valid_cluster = "0"

    # Create complete cluster per program
    folder_2_save_data = ""
    if (args.with_previous_common_clusters):
        folder_2_save_data += "COMMON"
    if (args.with_previous_individual_clusters):
        if ("COMMON" in folder_2_save_data):
            folder_2_save_data = "IND_" + folder_2_save_data
        else:
            folder_2_save_data = "IND"

    input_dir = os.path.join(args.root_input_dir, "Cfacial_clustering_results", args.dataset, folder_2_save_data,
                             args.assignation_method,
                             "assignation_" + args.assignation_type,
                             "clustering_th_" + str(args.qualityTh_valid_cluster),
                             "cluster_filtering_assignation")
    root_input_dir_parameters = os.path.join(args.root_input_dir, "Cfacial_clustering_results", args.dataset)

    output_dir = args.output_dir
    # Extract programs & participants
    extension = os.listdir(args.program_participants_folder)[0].split(".")[-1] \
        if ("." in os.listdir(args.program_participants_folder)[0]) else ''
    programs = sorted(os.listdir(input_dir)) if (
            args.program_name == None or args.program_name == 'None') else [args.program_name]

    participants_per_program = {program: list(set(['OTHER'] + loader.load_list_of_tertulianos(
                                                      args.program_participants_folder, program,
                                                      extension))) for program in programs}

    replace_underScore = True
    for program in programs:
        #Create objects per program:
        try:
            participants_per_program[program].remove("OTHER")
        except ValueError:
            pass
        input_path_mtcnn = os.path.join(args.root_path_MTCNN_results, program, "mtcnn_debug_sum")
        input_path_embs = os.path.join(args.root_path_MTCNN_results, program, "embeddings_sum")
        input_path_bbox = os.path.join(args.root_path_MTCNN_results, program, "boundingboxes_sum")
        input_path_models = os.path.join(args.root_input_dir, "Cfacial_clustering_results", args.dataset, "cluster_models",
                                         program)
        input_path_parameters = os.path.join(args.root_input_dir, "Cfacial_clustering_results", args.dataset,
                                             "cluster_graphics", "parameters", program)
        cleaner_obj = ClusterFeatures(input_dir, args.root_path_MTCNN_results, output_dir, program,
                                  participants_per_program[program],
                                  root_input_dir_parameters, input_path_models, input_path_parameters,
                                  args.first_param_name, args.second_param_name, args.quality_metric)
        for participant in sorted(participants_per_program[program]):
            if (os.path.exists(os.path.join(args.output_dir, participant.replace(" ", "_") + "_cluster_features_individual.csv"))):
                continue
            df_global_features, df_individual_clt_features = cleaner_obj.extract_clt_features(participant)
            df_global_features.to_csv(os.path.join(args.output_dir, participant.replace(" ", "_") + "_cluster_features_global.csv"), sep=";", index=False)
            df_individual_clt_features.to_csv(os.path.join(args.output_dir, participant.replace(" ", "_") + "_cluster_features_individual.csv"), sep=";",index=False)

ext: py

reason: import numpy


text:
#!/usr/bin/env python

"""Generators - Directed Graphs
----------------------------
"""

from nose.tools import *
from networkx_mod import *
from networkx_mod.generators.directed import *

class TestGeneratorsDirected():
    def test_smoke_test_random_graphs(self):
        G=gn_graph(100)
        G=gnr_graph(100,0.5)
        G=gnc_graph(100)
        G=scale_free_graph(100)

    def test_create_using_keyword_arguments(self):
        assert_raises(networkx_mod.exception.NetworkXError,
                      gn_graph, 100, create_using=Graph())
        assert_raises(networkx_mod.exception.NetworkXError,
                      gnr_graph, 100, 0.5, create_using=Graph())
        assert_raises(networkx_mod.exception.NetworkXError,
                      gnc_graph, 100, create_using=Graph())
        assert_raises(networkx_mod.exception.NetworkXError,
                      scale_free_graph, 100, create_using=Graph())
        G=gn_graph(100,seed=1)
        MG=gn_graph(100,create_using=MultiDiGraph(),seed=1)
        assert_equal(G.edges(), MG.edges())
        G=gnr_graph(100,0.5,seed=1)
        MG=gnr_graph(100,0.5,create_using=MultiDiGraph(),seed=1)
        assert_equal(G.edges(), MG.edges())
        G=gnc_graph(100,seed=1)
        MG=gnc_graph(100,create_using=MultiDiGraph(),seed=1)
        assert_equal(G.edges(), MG.edges())

reason: from networkx


text:
import numpy as np
from PIL import Image
import torch
from torch.utils.data import DataLoader, Dataset, ConcatDataset, Subset
from torchvision import datasets, transforms
from typing import Callable, Iterable, Tuple
from pathlib import Path

class NormalizeInverse(transforms.Normalize):
    """
    Undoes the normalization and returns the reconstructed images in the input domain.
    """

    def __init__(self, mean, std):
        mean = torch.as_tensor(mean)
        std = torch.as_tensor(std)
        std_inv = 1 / (std + 1e-7)
        mean_inv = -mean * std_inv
        super().__init__(mean=mean_inv, std=std_inv)

    def __call__(self, tensor):
        return super().__call__(tensor.clone())

CIFAR_PATH = Path("./data_cifar10")
CIFAR_TRANSFORM_NORMALIZE_MEAN = (0.4914, 0.4822, 0.4465)
CIFAR_TRANSFORM_NORMALIZE_STD = (0.2023, 0.1994, 0.2010)
CIFAR_TRANSFORM_NORMALIZE = transforms.Normalize(
    CIFAR_TRANSFORM_NORMALIZE_MEAN, CIFAR_TRANSFORM_NORMALIZE_STD
)
CIFAR_TRANSFORM_NORMALIZE_INV = NormalizeInverse(
    CIFAR_TRANSFORM_NORMALIZE_MEAN, CIFAR_TRANSFORM_NORMALIZE_STD
)
CIFAR_TRANSFORM_TRAIN = transforms.Compose(
    [
        transforms.RandomCrop(32, padding=4),
        transforms.RandomHorizontalFlip(),
        transforms.ToTensor(),
        CIFAR_TRANSFORM_NORMALIZE,
    ]
)
CIFAR_TRANSFORM_TRAIN_XY = lambda xy: (CIFAR_TRANSFORM_TRAIN(xy[0]), xy[1])

CIFAR_TRANSFORM_TEST = transforms.Compose(
    [
        transforms.ToTensor(),
        CIFAR_TRANSFORM_NORMALIZE,
    ]
)
CIFAR_TRANSFORM_TEST_XY = lambda xy: (CIFAR_TRANSFORM_TEST(xy[0]), xy[1])

class LabelSortedDataset(ConcatDataset):
    def __init__(self, dataset: Dataset):
        self.orig_dataset = dataset
        self.by_label = {}
        for i, (_, y) in enumerate(dataset):
            self.by_label.setdefault(y, []).append(i)

        self.n = len(self.by_label)
        assert set(self.by_label.keys()) == set(range(self.n))
        self.by_label = [Subset(dataset, self.by_label[i]) for i in range(self.n)]
        super().__init__(self.by_label)

    def subset(self, labels: Iterable[int]) -> ConcatDataset:
        if isinstance(labels, int):
            labels = [labels]
        return ConcatDataset([self.by_label[i] for i in labels])

class FilterDataset(Subset):
    def __init__(self, dataset: Dataset, *, label: int):
        indices = []
        for i, (_, y) in enumerate(dataset):
            if y == label:
                indices.append(i)
        super().__init__(dataset, indices)

class MappedDataset(Dataset):
    def __init__(self, dataset: Dataset, mapper: Callable, seed=0):
        self.dataset = dataset
        self.mapper = mapper
        self.seed = seed

    def __getitem__(self, i: int):
        if hasattr(self.mapper, 'seed'):
            self.mapper.seed(i + self.seed)
        return self.mapper(self.dataset[i])

    def __len__(self):
        return len(self.dataset)

class PoisonedDataset(Dataset):
    def __init__(
        self,
        dataset: Dataset,
        poisoner,
        *,
        label=None,
        indices=None,
        eps=500,
        seed=1,
        transform=None
    ):
        self.orig_dataset = dataset
        self.label = label
        if not indices and not eps:
            raise ValueError()

        if not indices:
            if label is not None:
                clean_inds = [i for i, (x, y) in enumerate(dataset) if y == label]
            else:
                clean_inds = range(len(dataset))

            rng = np.random.RandomState(seed)
            indices = rng.choice(clean_inds, eps, replace=False)

        self.indices = indices
        self.poison_dataset = MappedDataset(Subset(dataset, indices), poisoner, seed=seed)
        if transform:
            self.poison_dataset = MappedDataset(self.poison_dataset, transform)

        clean_indices = list(set(range(len(dataset))).difference(indices))
        self.clean_dataset = Subset(dataset, clean_indices)
        if transform:
            self.clean_dataset = MappedDataset(self.clean_dataset, transform)

        self.dataset = ConcatDataset([self.clean_dataset, self.poison_dataset])

    def __getitem__(self, i: int):
        return self.dataset[i]

    def __len__(self):
        return len(self.dataset)

class Poisoner(object):
    def poison(self, x: Image.Image) -> Image.Image:
        raise NotImplementedError()

    def __call__(self, x: Image.Image) -> Image.Image:
        return self.poison(x)

class PixelPoisoner(Poisoner):
    def __init__(
        self,
        *,
        method="pixel",
        pos: Tuple[int, int] = (11, 16),
        col: Tuple[int, int, int] = (101, 0, 25)
    ):
        self.method = method
        self.pos = pos
        self.col = col

    def poison(self, x: Image.Image) -> Image.Image:
        ret_x = x.copy()
        pos, col = self.pos, self.col

        if self.method == "pixel":
            ret_x.putpixel(pos, col)
        elif self.method == "pattern":
            ret_x.putpixel(pos, col)
            ret_x.putpixel((pos[0] - 1, pos[1] - 1), col)
            ret_x.putpixel((pos[0] - 1, pos[1] + 1), col)
            ret_x.putpixel((pos[0] + 1, pos[1] - 1), col)
            ret_x.putpixel((pos[0] + 1, pos[1] + 1), col)
        elif self.method == "ell":
            ret_x.putpixel(pos, col)
            ret_x.putpixel((pos[0] + 1, pos[1]), col)
            ret_x.putpixel((pos[0], pos[1] + 1), col)

        return ret_x

class StripePoisoner(Poisoner):
    def __init__(self, *, horizontal=True, strength=6, freq=16):
        self.horizontal = horizontal
        self.strength = strength
        self.freq = freq

    def poison(self, x: Image.Image) -> Image.Image:
        arr = np.asarray(x)
        (w, h, d) = arr.shape
        assert w == h  # have not tested w != h
        mask = np.full(
            (d, w, h), np.sin(np.linspace(0, self.freq * np.pi, h))
        ).swapaxes(0, 2)
        if self.horizontal:
            mask = mask.swapaxes(0, 1)
        mix = np.asarray(x) + self.strength * mask
        return Image.fromarray(np.uint8(mix.clip(0, 255)))

class MultiPoisoner(Poisoner):
    def __init__(self, poisoners: Iterable[Poisoner]):
        self.poisoners = poisoners

    def poison(self, x):
        for poisoner in self.poisoners:
            x = poisoner.poison(x)
        return x

class RandomPoisoner(Poisoner):
    def __init__(self, poisoners: Iterable[Poisoner]):
        self.poisoners = poisoners
        self.rng = np.random.RandomState()

    def poison(self, x):
        poisoner = self.rng.choice(self.poisoners)
        return poisoner.poison(x)

    def seed(self, i):
        self.rng.seed(i)

class LabelPoisoner(Poisoner):
    def __init__(self, poisoner: Poisoner, target_label: int):
        self.poisoner = poisoner
        self.target_label = target_label

    def poison(self, xy):
        x, _ = xy
        return self.poisoner(x), self.target_label

    def seed(self, i):
        if hasattr(self.poisoner, 'seed'):
            self.poisoner.seed(i)

def load_cifar_dataset(train=True):
    dataset = datasets.CIFAR10(root=str(CIFAR_PATH), train=train, download=True)
    return dataset

def make_dataloader(dataset: Dataset, batch_size, *, shuffle=True, drop_last=True):
    dataloader = DataLoader(
        dataset,
        batch_size=batch_size,
        shuffle=shuffle,
        num_workers=4,
        pin_memory=True,
        drop_last=drop_last,
    )
    return dataloader

def load_cifar_train(batch_size=32):
    path = "./data_cifar10"
    kwargs = {"num_workers": 4, "pin_memory": True, "drop_last": True}
    transform_train = transforms.Compose(
        [
            transforms.RandomCrop(32, padding=4),
            transforms.RandomHorizontalFlip(),
            transforms.ToTensor(),
            transforms.Normalize((0.4914, 0.4822, 0.4465), (0.2023, 0.1994, 0.2010)),
        ]
    )
    trainset = datasets.CIFAR10(
        root=path, train=True, download=True, transform=transform_train
    )
    trainloader = DataLoader(trainset, batch_size=batch_size, shuffle=True, **kwargs)
    return trainloader

def load_cifar_test(batch_size=32):
    path = "./data_cifar10"
    kwargs = {"num_workers": 4, "pin_memory": True, "drop_last": True}
    transform_test = transforms.Compose(
        [
            transforms.ToTensor(),
            transforms.Normalize((0.4914, 0.4822, 0.4465), (0.2023, 0.1994, 0.2010)),
        ]
    )
    testset = datasets.CIFAR10(
        root=path, train=False, download=True, transform=transform_test
    )
    testloader = DataLoader(testset, batch_size=batch_size, shuffle=False, **kwargs)
    return testloader

max_stars_repo_name: TDteach/spectre-defense

max_issues_repo_name: TDteach/spectre-defense

max_forks_repo_name: TDteach/spectre-defense


text:
"""Objective functions for hyperparameter optimizations without training."""
import numpy as np
import torch
from torch import Tensor
from torch.utils.data import DataLoader

from mala.common.printout import printout
from mala.common.parameters import Parameters
from mala.datahandling.data_handler import DataHandler
from mala.network.network import Network
from mala.network.objective_base import ObjectiveBase

class ObjectiveNoTraining(ObjectiveBase):
    """
    Represents the objective function using no NN training.

    The objective value is calculated using the Jacobian.

    Parameters
    ----------
    search_parameters : mala.common.parameters.Parameters
        Parameters used to create this objective.

    data_handler : mala.datahandling.data_handler.DataHandler
        datahandler to be used during the hyperparameter optimization.

    trial_type : string
        Format of hyperparameters used in this objective. Supported
        choices are:

            - optuna
            - oat
    """

    def __init__(self, search_parameters: Parameters, data_handler:
                 DataHandler, trial_type):
        super(ObjectiveNoTraining, self).__init__(search_parameters,
                                                  data_handler)
        self.trial_type = trial_type

    def __call__(self, trial):
        """
        Get objective value for a trial (=set of hyperparameters).

        Parameters
        ----------
        trial
            A trial is a set of hyperparameters; can be an optuna based
            trial or simply a OAT compatible list.
        """
        # Parse the parameters using the base class.
        super(ObjectiveNoTraining, self).parse_trial(trial)

        # Build the network.
        surrogate_losses = []
        for i in range(0, self.params.hyperparameters.
                       number_training_per_trial):
            net = Network(self.params)
            device = "cuda" if self.params.use_gpu else "cpu"

            # Load the batchesand get the jacobian.
            do_shuffle = self.params.running.use_shuffling_for_samplers
            if self.data_handler.parameters.use_lazy_loading or \
                    self.params.use_horovod:
                do_shuffle = False
            if self.params.running.use_shuffling_for_samplers:
                self.data_handler.mix_datasets()
            loader = DataLoader(self.data_handler.training_data_set,
                                batch_size=self.params.running.mini_batch_size,
                                shuffle=do_shuffle)
            jac = ObjectiveNoTraining.__get_batch_jacobian(net, loader, device)

            # Loss = - score!
            surrogate_loss = float('inf')
            try:
                surrogate_loss = - ObjectiveNoTraining.__calc_score(jac)
                surrogate_loss = surrogate_loss.cpu().detach().numpy().astype(
                    np.float64)
            except RuntimeError:
                printout("Got a NaN, ignoring sample.")
            surrogate_losses.append(surrogate_loss)

        if self.params.hyperparameters.number_training_per_trial > 1:
            printout("Losses from multiple runs are: ")
            printout(surrogate_losses)

        if self.params.hyperparameters.trial_ensemble_evaluation == "mean":
            return np.mean(surrogate_losses)

        elif self.params.hyperparameters.trial_ensemble_evaluation == \
                "mean_std":
            mean = np.mean(surrogate_losses)

            # Cannot calculate the standar deviation of a bunch of infinities.
            if np.isinf(mean):
                return mean
            else:
                return np.mean(surrogate_losses) + \
                       np.std(surrogate_losses)
        else:
            raise Exception("No way to estimate the trial metric from ensemble"
                            " training provided.")

    @staticmethod
    def __get_batch_jacobian(net: Network, loader: DataLoader, device) \
            -> Tensor:
        """Calculate the jacobian of the batch."""
        x: Tensor
        (x, _) = next(iter(loader))
        x = x.to(device)
        net.zero_grad()
        x.requires_grad_(True)
        y: Tensor = net(x)
        y.backward(torch.ones_like(y))
        jacobian = x.grad.detach()
        return jacobian

    @staticmethod
    def __corrcoef(x):
        """
        Calculate the correlation coefficient of an array.

        Pytorch implementation of np.corrcoef.
        """
        mean_x = torch.mean(x, 1, True)
        xm = x.sub(mean_x.expand_as(x))
        c = xm.mm(xm.t())
        c = c / (x.size(1) - 1)

        # normalize covariance matrix
        d = torch.diag(c)
        stddev = torch.pow(d, 0.5)
        c = c.div(stddev.expand_as(c))
        c = c.div(stddev.expand_as(c).t())

        # clamp between -1 and 1
        # probably not necessary but numpy does it
        c = torch.clamp(c, -1.0, 1.0)

        return c

    @staticmethod
    def __calc_score(jacobian: Tensor):
        """Calculate the score for a Jacobian."""
        correlations = ObjectiveNoTraining.__corrcoef(jacobian)
        eigen_values, _ = torch.eig(correlations)
        # Only consider the real valued part, imaginary part is rounding
        # artefact
        eigen_values = eigen_values.type(torch.float32)
        # Needed for numerical stability. 1e-4 instead of 1e-5 in reference
        # as the torchs eigenvalue solver on GPU
        # seems to have bigger rounding errors than numpy, resulting in
        # slightly larger negative Eigenvalues
        k = 1e-4
        v = -torch.sum(torch.log(eigen_values + k) + 1. / (eigen_values+k))
        return v

max_stars_repo_name: DanielKotik/mala

max_issues_repo_name: DanielKotik/mala

max_forks_repo_name: DanielKotik/mala


text:
# -*- coding: utf-8 -*-
"""
Created on Wed Jul  8 14:35:31 2020

@author: sopmathieu

This file contains functions to apply the monitoring procedure (composed of
a CUSUM chart and support vector classifier and regressor) on a time-series.
The script also contains functions to graphically display the results of
this monitoring.

"""

import numpy as np
import matplotlib.pyplot as plt

from cusvm import svm_training as svm

def firstNonNan(x):
  """
  Finds first non NaN value in an array.

  Parameters
  ----------
  x : 1D-array
     Initial array to retrieve its first non NaN value.

  Returns
  ------
  ind : int
      Index of first non NaN value of the array.

  """
  for ind in range(len(x)):
    if not np.isnan(x[ind]):
      return ind

def alerts_info(data, L_plus, delta, wdw_length, clf, reg, L_minus=None,
                     k=None, cut=None, verbose=True, wdw_shift=0):
    """
    Applies the two-sided CUSUM chart on a series.

    This function returns the size and the form ('jumps', 'drifts', 'oscillating
    shifts') of the shifts after each alert.

    Parameters
    ----------
    data : 1D-array
        A single series of standardized observations to be monitored.
    L_plus : float
        Value for the positive control limit.
    delta : float >= 0
        The target shift size.
    wdw_length : int > 0
        The length of the input vector.
    clf : support vector classification model
       The trained classifier.
    reg : support vector regression model
        The trained regressor.
    L_minus :  float, optional
        Value for the negative control limit. Default is None.
        When None, L_minus = - L_plus.
    k : float, optional
        The allowance parameter. The default is None.
        When None, k = delta/2 (optimal formula for iid normal data).
    cut : float, optional
        Upper value for the chart statistics. Values of the positive (resp.
        negative) chart statistics are constrained to be equal to or
        lower than 'cut' (resp. equal to or superior than -'cut').
        When None, cut is equal to '2L_plus'. The default is None.
    Verbose : bool, optional
        Flag to print the percentage of alert in the series. Default is True.
    wdw_shift : int, optional
        Shift that is applied to the predictions of the SVMs, to be better
        aligned with the actual deviations.
        The default is 0.

    Returns
    -------
    form_plus : 1D-array
        Predicted shift forms after positive alerts.
        When no alerts is detected, the shift forms are set to NaNs.
    form_minus : 1D-array
        Predicted shift forms after negative alerts.
        When no alerts is detected, the shift forms are set to NaNs.
    size_plus : 1D-array
        Predicted shift sizes after positive alerts.
        When no alerts is detected, the shift sizes are set to NaNs.
    size_minus : 1D-array
        Predicted shift sizes after negative alerts.
        When no alerts is detected, the shift sizes are set to NaNs.
    C_plus :
        Values of the positive chart statistic.
    C_minus :
        Values of the negative chart statistic.

    """
    assert np.ndim(data) == 1, "Input data must be a 1D array (one series)"
    n_obs = len(data)
    if L_minus is None:
        L_minus = -L_plus
    if k is None:
        k = abs(delta)/2
    if cut is None:
        cut = L_plus * 2
    wdw_shift = int(wdw_shift)
    assert wdw_shift < wdw_length, "wdw_shift should be inferior to wdw_length"

    length = len(data[~np.isnan(data)])

    input_minus = np.zeros((n_obs, wdw_length)); input_plus = np.zeros((n_obs, wdw_length))
    input_minus[:] = np.nan; input_plus[:] = np.nan
    flag_plus = np.zeros((n_obs)); flag_minus = np.zeros((n_obs))
    C_plus = np.zeros((n_obs)); C_minus = np.zeros((n_obs))

    for i in range(wdw_length, n_obs):

        ## CUSUM monitoring
        C_plus[i] =  min(cut, max(0, C_plus[i-1] + data[i] - k))#avoid cusum "explosion"
        if C_plus[i] > L_plus: #alert
            flag_plus[i] = 1
            input_plus[i,:] = data[i+1-wdw_length:i+1]

        C_minus[i] = max(-cut, min(0, C_minus[i-1] + data[i] + k))
        if C_minus[i] < L_minus: #alert
            flag_minus[i] = 1
            input_minus[i,:] = data[i+1-wdw_length:i+1]

    ## compute percentage of alerts
    oc_p = np.nonzero(flag_plus)
    oc_m = np.nonzero(flag_minus)
    #if alert both for pos and neg limits, count for only one alert
    oc_both = len(set(np.concatenate((oc_p[0], oc_m[0]))))
    #OC_perc = oc_both*100/n_obs #total period
    OC_perc = oc_both*100/length #observing period

    if verbose :
        print("Percentage of alerts: %0.2f" %OC_perc)

    ## interpolate NaNs in input vectors
    input_minus_valid, ind_minus = svm.fill_nan(input_minus)
    input_plus_valid, ind_plus = svm.fill_nan(input_plus)

    ## shift the indexes of the SVM predictions
    ## otherwise the deviation predicted at time t is based on
    ## the entire wdw_length
    ind_minus = ind_minus - wdw_shift
    ind_plus = ind_plus - wdw_shift

    ##apply classifier and regressor on (filled-up) input vectors
    size_minus = np.zeros((n_obs)); size_plus = np.zeros((n_obs))
    size_minus[:] = np.nan; size_plus[:] = np.nan
    form_minus = np.zeros((n_obs)); form_plus = np.zeros((n_obs))
    form_minus[:] = np.nan; form_plus[:] = np.nan

    if len(ind_minus)>0: #at least one value
        size_minus[ind_minus] = reg.predict(input_minus_valid)
        form_minus[ind_minus] = clf.predict(input_minus_valid)
    if len(ind_plus)>0:
        size_plus[ind_plus] = reg.predict(input_plus_valid)
        form_plus[ind_plus] = clf.predict(input_plus_valid)

    return (form_plus, form_minus, size_plus, size_minus, C_plus, C_minus)

#==================================================================
#==================================================================

def plot_4panels(data, unstn_data, L_plus, time, form_plus, form_minus, size_plus,
                      size_minus, C_plus, C_minus, name, L_minus=None,
                      time_start=None, time_stop=None, x_ticks=None,
                      labels_ticks=None, fig_size=[7,10]):
    """
    Plots the main results of the monitoring.

    This function plots different quantities that resume the main features
    of the monitoring into four panels.
    The first panel represents the unstandadised data,
    the second shows the standardized data,
    the third panel displays the CUSUM statistics applied to the standardised data
    in square-root scale
    and the last panel shows the predicted shift sizes and forms.

    Parameters
    ---------
    data : 1D-array
        A (single) series of standardized observations to be monitored.
    unstn_data : 1D-array
        The original series (without standardization).
    L_plus : float
        Value for the positive control limit.
    time : 1D-array
        An array with the time of the observations (same length as data).
    form_plus : 1D-array
        Predicted shift forms after positive alerts (same length as data).
        When no alerts is detected, the shift forms should be set to NaNs.
    form_minus : 1D-array
        Predicted shift forms after negative alerts (same length as data).
        When no alerts is detected, the shift forms should be set to NaNs.
    size_plus : 1D-array
        Predicted shift sizes after positive alerts (same length as data).
        When no alerts is detected, the shift sizes should be set to NaNs.
    size_minus : 1D-array
        Predicted shift sizes after negative alerts (same length as data).
        When no alerts is detected, the shift sizes should be set to NaNs.
    C_plus : 1D-array
        Values of the positivie chart statistic (same length as data).
    C_minus : 1D-array
        Values of the negative chart statistic (same length as data).
    name : str
       The index of the series.
    L_minus :  float, optional
        Value for the negative control limit. Default is None.
        When None, L_minus = - L_plus.
    time_start : int, optional
       Starting time value of the plot, in years. The plot shows the period
       [time_start, time_stop]. Default is None.
       When None, the first non-Nan observation in the series is used.
    time_stop : int, optional
      Last time value of the plot, in years. The plot shows the period
      [time_start, time_stop]. Default is None.
      When None, the series is represented until the last observation.
    x_ticks : 1D-array, optional
        Locations of the ticks for the x-axis of the plot.
        The default is None.
    labels_ticks : 1D-array, optional
        Labels of the ticks for the x-axis of the plot.
        The labels can only be passed if the locations are passed as well.
        The default is None.
    fig_size : list, optional
        Width and heigth of the resulting figure.
        The default is [7,10].

    Returns
    ------
    fig : a matplotlib figure
        The figure (with four panels).

    """

    assert np.ndim(data) == 1, "Input data must be a 1D array (one series)"
    n_obs = len(data)
    if time_start is None:
        start = firstNonNan(data)
    else:
        try:
             time[-1] >= time_start
        except TypeError as e:
            e.args += ("the values of time should be of the same type as time_start",)
            raise
        start = [i for i in range(len(time)) if time[i] >= time_start][0]
    if time_stop is None:
        stop = len(time) - 1
    else:
        try:
             time[-1] >= time_stop
        except TypeError as e:
            e.args += ("the values of time should be of the same type as time_stop",)
            raise
        stop = [i for i in range(len(time)) if time[i] >= time_stop][0]
    assert stop > start, "time_stop should be stricly superior to time_start!"

    if L_minus is None:
        L_minus = -L_plus

    #colors
    colorInd_plus = np.where(~np.isnan(form_plus))[0][:]
    colorInd_minus = np.where(~np.isnan(form_minus))[0][:]
    color_graph_p = np.ones((n_obs))*3; color_graph_m = np.ones((n_obs))*3
    color_graph_m[colorInd_minus] = form_minus[colorInd_minus]
    color_graph_p[colorInd_plus] = form_plus[colorInd_plus]

    #jumps
    jump_p = np.zeros((n_obs)); jump_p[:]=np.nan
    jump_m = np.zeros((n_obs)); jump_m[:]=np.nan
    jump_m[np.where(color_graph_m == 0)[0]] = size_minus[np.where(color_graph_m == 0)[0]]
    jump_p[np.where(color_graph_p == 0)[0]] = size_plus[np.where(color_graph_p == 0)[0]]

    #trends
    trend_p = np.zeros((n_obs)); trend_p[:] = np.nan
    trend_m = np.zeros((n_obs)); trend_m[:] = np.nan
    trend_m[np.where(color_graph_m == 1)[0]] = size_minus[np.where(color_graph_m == 1)[0]]
    trend_p[np.where(color_graph_p == 1)[0]] = size_plus[np.where(color_graph_p == 1)[0]]

    #oscillating shifts
    oscill_p = np.zeros((n_obs)); oscill_p[:] = np.nan
    oscill_m = np.zeros((n_obs)); oscill_m[:] = np.nan
    oscill_m[np.where(color_graph_m == 2)[0]] = size_minus[np.where(color_graph_m == 2)[0]]
    oscill_p[np.where(color_graph_p == 2)[0]] = size_plus[np.where(color_graph_p == 2)[0]]

    plt.rcParams['font.size'] = 14

    data_per = data[start:stop]
    max_val = max(data_per[~np.isnan(data_per)])
    min_val = min(data_per[~np.isnan(data_per)])
    y_max = 1.2*max_val ; y_min = 1.2*min_val

    fig = plt.figure(figsize=(fig_size[0],fig_size[1]))
    f1 = fig.add_subplot(4, 1, 1)
    plt.title("Monitoring in %s" %name)
    plt.ylabel('Unstand. residuals')
    plt.plot(time[start:stop], unstn_data[start:stop])
    plt.plot([time[start], time[stop]], [1, 1], 'k-', lw=2)
    f1.axes.get_xaxis().set_ticklabels([])
    f1.set_xlim([time[start], time[stop]])
    if x_ticks is not None:
        plt.xticks(x_ticks)
    f1.axes.xaxis.grid(True, linewidth=0.15)

    f2 = fig.add_subplot(4, 1, 2)
    plt.ylabel('Stand. residuals')
    plt.plot(time[start:stop], data[start:stop])
    plt.plot([time[start], time[stop]], [0, 0], 'k-', lw=2)
    f2.set_ylim([y_min, y_max]); f2.set_xlim([time[start], time[stop]])
    if x_ticks is not None:
        plt.xticks(x_ticks)
    f2.axes.get_xaxis().set_ticklabels([])
    f2.axes.xaxis.grid(True, linewidth=0.15)

    f3 = fig.add_subplot(4, 1, 3)
    plt.ylabel('CUSUM stats')
    plt.plot(time[start:stop], np.sqrt(C_plus[start:stop]), c='tab:red', label='$C^+$')
    plt.plot(time[start:stop], -np.sqrt(-C_minus[start:stop]), '--', c='tab:brown', label='$C^-$')
    #upper control limit value (horizontal line)
    plt.plot([time[start], time[stop]], [np.sqrt(L_plus), np.sqrt(L_plus)], 'k-', lw=2.2)
    #lower control limit value (horizontal line)
    plt.plot([time[start], time[stop]], [-np.sqrt(abs(L_minus)), -np.sqrt(abs(L_minus))], 'k-', lw=2.2)
    plt.legend(loc='upper right', ncol=2, fontsize=10)
    f3.set_ylim([-np.sqrt(2*L_plus)*2, np.sqrt(2*L_plus)*2])
    f3.set_xlim([time[start], time[stop]])
    if x_ticks is not None:
        plt.xticks(x_ticks)
    f3.axes.get_xaxis().set_ticklabels([])
    f3.axes.xaxis.grid(True, linewidth=0.15)

    # sizes = np.array([size_plus, size_minus])
    # sizes_per = sizes[:,start:stop]
    # if len(sizes_per[~np.isnan(sizes_per)]) > 0:
    #     max_val = max(sizes_per[~np.isnan(sizes_per)])
    #     min_val = min(sizes_per[~np.isnan(sizes_per)])
    #     y_max = 1.5*max_val ; y_min = 1.5*min_val

    f4 = fig.add_subplot(4, 1, 4)
    plt.ylabel('Deviations')
    plt.plot(time[start:stop], data[start:stop], color='#7f7f7f', lw=0.2)
    plt.plot(time[start:stop], jump_m[start:stop], '--', c='tab:purple', label='jumps')
    plt.plot(time[start:stop], jump_p[start:stop], '--', c='tab:purple')
    plt.plot(time[start:stop], trend_m[start:stop],  c='tab:green', label='trends')
    plt.plot(time[start:stop], trend_p[start:stop],  c='tab:green')
    plt.plot(time[start:stop], oscill_m[start:stop], ':', c='tab:orange', label='oscill')
    plt.plot(time[start:stop], oscill_p[start:stop], ':', c='tab:orange')
    plt.plot([time[start], time[stop]], [0, 0], 'k-', lw=2)
    plt.legend(loc='lower right', ncol=3, fontsize=10)
    f4.set_ylim([y_min, y_max]); f4.set_xlim([time[start], time[stop]])
    if x_ticks is not None:
        if labels_ticks is not None:
            plt.xticks(x_ticks, labels_ticks)
        else:
            plt.xticks(x_ticks)
    plt.xlabel('time')
    plt.tick_params(axis='x', which='major')
    f4.axes.xaxis.grid(True, linewidth=0.15)

    #plt.tight_layout()
    plt.show()
    return fig

#==================================================================
### other plots
#==================================================================

def plot_1panel(data, time, form_plus, form_minus, size_plus,
             size_minus, name, time_start=None, time_stop=None,
             x_ticks=None, labels_ticks=None, fig_size=[12,4]):
    """
    Plots the predictions (forms and sizes) of the support
    vector machines (SVM).

    Parameters
    ---------
    data : 1D-array
        A (single) series of standardized observations to be monitored.
    time : 1D-array
        An array with the time of the observations (same length as data).
    form_plus : 1D-array
        Predicted shift forms after positive alerts (same length as data).
        When no alerts is detected, the shift forms should be set to NaNs.
    form_minus : 1D-array
        Predicted shift forms after negative alerts (same length as data).
        When no alerts is detected, the shift forms should be set to NaNs.
    size_plus : 1D-array
        Predicted shift sizes after positive alerts (same length as data).
        When no alerts is detected, the shift sizes should be set to NaNs.
    size_minus : 1D-array
        Predicted shift sizes after negative alerts (same length as data).
        When no alerts is detected, the shift sizes should be set to NaNs.
    name : str
       The index or the name of the series.
    time_start : int, optional
       Starting time value of the plot, in years. The plot shows the period
       [time_start, time_stop]. Default is None.
       When None, the first non-Nan observation in the series is used.
    time_stop : int, optional
      Last time value of the plot, in years. The plot shows the period
      [time_start, time_stop]. Default is None.
      When None, the series is represented until the last observation.
    x_ticks : 1D-array, optional
        Locations of the ticks for the x-axis of the plot.
        The default is None.
    labels_ticks : 1D-array, optional
        Labels of the ticks for the x-axis of the plot.
        The labels can only be passed if the locations are passed as well.
        The default is None.
    fig_size : list, optional
        Width and heigth of the resulting figure.
        The default is [12,4].

    Returns
    ------
    fig : a matplotlib figure
        The figure with the SVM predictions.

    """

    assert np.ndim(data) == 1, "Input data must be a 1D array (one series)"
    n_obs = len(data)
    if time_start is None:
        start = firstNonNan(data)
    else:
        try:
             time[-1] >= time_start
        except TypeError as e:
            e.args += ("the values of time should be of the same type as time_start",)
            raise
        start = [i for i in range(len(time)) if time[i] >= time_start][0]
    if time_stop is None:
        stop = len(time) - 1
    else:
        try:
             time[-1] >= time_stop
        except TypeError as e:
            e.args += ("the values of time should be of the same type as time_stop",)
            raise
        stop = [i for i in range(len(time)) if time[i] >= time_stop][0]

    assert stop > start, "time_stop should be stricly superior to time_start!"

    #colors
    colorInd_plus = np.where(~np.isnan(form_plus))[0][:]
    colorInd_minus = np.where(~np.isnan(form_minus))[0][:]
    color_graph_p = np.ones((n_obs))*3; color_graph_m = np.ones((n_obs))*3
    color_graph_m[colorInd_minus] = form_minus[colorInd_minus]
    color_graph_p[colorInd_plus] = form_plus[colorInd_plus]

    #jumps
    jump_p = np.zeros((n_obs)); jump_p[:]=np.nan
    jump_m = np.zeros((n_obs)); jump_m[:]=np.nan
    jump_m[np.where(color_graph_m == 0)[0]] = size_minus[np.where(color_graph_m == 0)[0]]
    jump_p[np.where(color_graph_p == 0)[0]] = size_plus[np.where(color_graph_p == 0)[0]]

    #trends
    trend_p = np.zeros((n_obs)); trend_p[:] = np.nan
    trend_m = np.zeros((n_obs)); trend_m[:] = np.nan
    trend_m[np.where(color_graph_m == 1)[0]] = size_minus[np.where(color_graph_m == 1)[0]]
    trend_p[np.where(color_graph_p == 1)[0]] = size_plus[np.where(color_graph_p == 1)[0]]

    #oscillating shifts
    oscill_p = np.zeros((n_obs)); oscill_p[:] = np.nan
    oscill_m = np.zeros((n_obs)); oscill_m[:] = np.nan
    oscill_m[np.where(color_graph_m == 2)[0]] = size_minus[np.where(color_graph_m == 2)[0]]
    oscill_p[np.where(color_graph_p == 2)[0]] = size_plus[np.where(color_graph_p == 2)[0]]

    plt.rcParams['font.size'] = 14

    sizes = np.array([size_plus, size_minus])
    sizes_per = sizes[:,start:stop]
    if len(sizes_per[~np.isnan(sizes_per)]) > 0:
        max_val = max(sizes_per[~np.isnan(sizes_per)])
        min_val = min(sizes_per[~np.isnan(sizes_per)])
        y_max = 1.5*max_val ; y_min = 1.5*min_val
    else:
        data_per = data[start:stop]
        max_val = max(data_per[~np.isnan(data_per)])
        min_val = min(data_per[~np.isnan(data_per)])
        y_max = 1.2*max_val ; y_min = 1.2*min_val

    fig = plt.figure(figsize=(fig_size[0],fig_size[1]))
    f1 = fig.add_subplot(1, 1, 1)
    plt.ylabel('Deviations')
    plt.title("Monitoring in %s" %name)
    plt.plot(time[start:stop], data[start:stop], color='#7f7f7f', lw=0.2)
    plt.plot(time[start:stop], jump_m[start:stop], '--', c='tab:purple', label='jumps')
    plt.plot(time[start:stop], jump_p[start:stop], '--', c='tab:purple')
    plt.plot(time[start:stop], trend_m[start:stop],  c='tab:green', label='trends')
    plt.plot(time[start:stop], trend_p[start:stop],  c='tab:green')
    plt.plot(time[start:stop], oscill_m[start:stop], ':', c='tab:orange', label='oscill')
    plt.plot(time[start:stop], oscill_p[start:stop], ':', c='tab:orange')
    plt.plot([time[start], time[stop]], [0, 0], 'k-', lw=2)
    plt.legend(loc='lower center', ncol=3)
    f1.set_ylim([y_min, y_max]); f1.set_xlim([time[start], time[stop]])
    if x_ticks is not None:
        if labels_ticks is not None:
            plt.xticks(x_ticks, labels_ticks)
        else:
            plt.xticks(x_ticks)
    plt.xlabel('time')
    plt.tick_params(axis='x', which='major')
    f1.axes.xaxis.grid(True, linewidth=0.15)

    #plt.tight_layout()
    plt.show()
    return fig

def plot_3panels(data, L_plus, time, form_plus, form_minus, size_plus,
                      size_minus, C_plus, C_minus, name, L_minus=None,
                      time_start=None, time_stop=None, x_ticks=None,
                      labels_ticks=None, fig_size=[7,7]):
    """
    Plots the main results of the monitoring.

    This function plots different quantities that resume the main features
    of the monitoring into three panels.
    The first panel represents the standardized data
    the second panel displays the CUSUM statistics applied to the standardised
    data in square-root scale
    and the last panel shows the predicted shift sizes and forms.

    Parameters
    ---------
    data : 1D-array
        A (single) series of standardized observations to be monitored.
    L_plus : float
        Value for the positive control limit.
    time : 1D-array
        An array with the time of the observations (same length as data).
    form_plus : 1D-array
        Predicted shift forms after positive alerts (same length as data).
        When no alerts is detected, the shift forms should be set to NaNs.
    form_minus : 1D-array
        Predicted shift forms after negative alerts (same length as data).
        When no alerts is detected, the shift forms should be set to NaNs.
    size_plus : 1D-array
        Predicted shift sizes after positive alerts (same length as data).
        When no alerts is detected, the shift sizes should be set to NaNs.
    size_minus : 1D-array
        Predicted shift sizes after negative alerts (same length as data).
        When no alerts is detected, the shift sizes should be set to NaNs.
    C_plus : 1D-array
        Values of the positivie chart statistic (same length as data).
    C_minus : 1D-array
        Values of the negative chart statistic (same length as data).
    name : str
       The index or the name of the series.
    L_minus :  float, optional
        Value for the negative control limit. Default is None.
        When None, L_minus = - L_plus.
    time_start : float or datetime object, optional
       Starting time value of the plot. Should be
       The plot shows the period
       [time_start, time_stop]. Default is None.
       When None, the first non-Nan observation in the series is used.
    time_stop : float or datetime object, optional
      Last time value of the plot. The plot shows the period
      [time_start, time_stop]. Default is None.
      When None, the series is represented until the last observation.
    x_ticks : 1D-array, optional
        Locations of the ticks for the x-axis of the plot.
        The default is None.
    labels_ticks : 1D-array, optional
        Labels of the ticks for the x-axis of the plot.
        The labels can only be passed if the locations are passed as well.
        The default is None.
    fig_size : list, optional
        Width and heigth of the resulting figure.
        The default is [7,7].

    Returns
    ------
    fig : a matplotlib figure
        The figure (with four panels).

    """

    assert np.ndim(data) == 1, "Input data must be a 1D array (one series)"
    n_obs = len(data)

    if time_start is None:
        start = firstNonNan(data)
    else:
        try:
             time[-1] >= time_start
        except TypeError as e:
            e.args += ("the values of time should be of the same type as time_start",)
            raise
        start = [i for i in range(len(time)) if time[i] >= time_start][0]
    if time_stop is None:
        stop = len(time) - 1
    else:
        try:
             time[-1] >= time_stop
        except TypeError as e:
            e.args += ("the values of time should be of the same type as time_stop",)
            raise
        stop = [i for i in range(len(time)) if time[i] >= time_stop][0]

    assert stop > start, "time_stop should be stricly superior to time_start!"

    if L_minus is None:
        L_minus = -L_plus

    #colors
    colorInd_plus = np.where(~np.isnan(form_plus))[0][:]
    colorInd_minus = np.where(~np.isnan(form_minus))[0][:]
    color_graph_p = np.ones((n_obs))*3; color_graph_m = np.ones((n_obs))*3
    color_graph_m[colorInd_minus] = form_minus[colorInd_minus]
    color_graph_p[colorInd_plus] = form_plus[colorInd_plus]

    #jumps
    jump_p = np.zeros((n_obs)); jump_p[:]=np.nan
    jump_m = np.zeros((n_obs)); jump_m[:]=np.nan
    jump_m[np.where(color_graph_m == 0)[0]] = size_minus[np.where(color_graph_m == 0)[0]]
    jump_p[np.where(color_graph_p == 0)[0]] = size_plus[np.where(color_graph_p == 0)[0]]

    #trends
    trend_p = np.zeros((n_obs)); trend_p[:] = np.nan
    trend_m = np.zeros((n_obs)); trend_m[:] = np.nan
    trend_m[np.where(color_graph_m == 1)[0]] = size_minus[np.where(color_graph_m == 1)[0]]
    trend_p[np.where(color_graph_p == 1)[0]] = size_plus[np.where(color_graph_p == 1)[0]]

    #oscillating shifts
    oscill_p = np.zeros((n_obs)); oscill_p[:] = np.nan
    oscill_m = np.zeros((n_obs)); oscill_m[:] = np.nan
    oscill_m[np.where(color_graph_m == 2)[0]] = size_minus[np.where(color_graph_m == 2)[0]]
    oscill_p[np.where(color_graph_p == 2)[0]] = size_plus[np.where(color_graph_p == 2)[0]]

    plt.rcParams['font.size'] = 14

    data_per = data[start:stop]
    max_val = max(data_per[~np.isnan(data_per)])
    min_val = min(data_per[~np.isnan(data_per)])
    y_max = 1.2*max_val ; y_min = 1.2*min_val

    fig = plt.figure(figsize=(fig_size[0],fig_size[1]))
    f1 = fig.add_subplot(3, 1, 1)
    plt.title("Monitoring in %s" %name)
    plt.ylabel('Stand. residuals')
    plt.plot(time[start:stop], data[start:stop])
    plt.plot([time[start], time[stop]], [0, 0], 'k-', lw=2)
    f1.set_ylim([y_min, y_max]); f1.set_xlim([time[start], time[stop]])
    if x_ticks is not None:
        plt.xticks(x_ticks)
    f1.axes.get_xaxis().set_ticklabels([])
    f1.axes.xaxis.grid(True, linewidth=0.1, linestyle='-', color='#7f7f7f')

    f2 = fig.add_subplot(3, 1, 2)
    plt.ylabel('CUSUM stats')
    plt.plot(time[start:stop], np.sqrt(C_plus[start:stop]), c='tab:red', label='$C^+$')
    plt.plot(time[start:stop], -np.sqrt(-C_minus[start:stop]), '--', c='tab:brown', label='$C^-$')
    #upper control limit value (horizontal line)
    plt.plot([time[start], time[stop]], [np.sqrt(L_plus), np.sqrt(L_plus)], 'k-', lw=2.2)
    #lower control limit value (horizontal line)
    plt.plot([time[start], time[stop]], [-np.sqrt(abs(L_minus)), -np.sqrt(abs(L_minus))], 'k-', lw=2.2)
    plt.legend(loc='upper right', ncol=2, fontsize=10)
    f2.set_ylim([-np.sqrt(2*L_plus)*2, np.sqrt(2*L_plus)*2])
    f2.set_xlim([time[start], time[stop]])
    if x_ticks is not None:
        plt.xticks(x_ticks)
    f2.axes.get_xaxis().set_ticklabels([])
    f2.axes.xaxis.grid(True, linewidth=0.1, linestyle='-', color='#7f7f7f')

    # sizes = np.array([size_plus, size_minus])
    # sizes_per = sizes[:,start:stop]
    # if len(sizes_per[~np.isnan(sizes_per)]) > 0:
    #     max_val = max(sizes_per[~np.isnan(sizes_per)])
    #     min_val = min(sizes_per[~np.isnan(sizes_per)])
    #     y_max = 1.5*max_val ; y_min = 1.5*min_val

    f3 = fig.add_subplot(3, 1, 3)
    plt.ylabel('Deviations')
    plt.plot(time[start:stop], data[start:stop], color='#7f7f7f', lw=0.2)
    plt.plot(time[start:stop], jump_m[start:stop], '--', c='tab:purple', label='jumps')
    plt.plot(time[start:stop], jump_p[start:stop], '--', c='tab:purple')
    plt.plot(time[start:stop], trend_m[start:stop],  c='tab:green', label='trends')
    plt.plot(time[start:stop], trend_p[start:stop],  c='tab:green')
    plt.plot(time[start:stop], oscill_m[start:stop], ':', c='tab:orange', label='oscill')
    plt.plot(time[start:stop], oscill_p[start:stop], ':', c='tab:orange')
    plt.plot([time[start], time[stop]], [0, 0], 'k-', lw=2)
    plt.legend(loc='lower right', ncol=3, fontsize=10)
    f3.set_ylim([y_min, y_max])
    f3.set_xlim([time[start], time[stop]])
    if x_ticks is not None:
        if labels_ticks is not None:
            plt.xticks(x_ticks, labels_ticks)
        else:
            plt.xticks(x_ticks)
    plt.xlabel('time')
    plt.tick_params(axis='x', which='major')
    f3.axes.xaxis.grid(True, linewidth=0.1, linestyle='-', color='#7f7f7f')

    #plt.tight_layout()
    plt.show()

    return fig

max_stars_repo_name: sophiano/cusvm

max_issues_repo_name: sophiano/cusvm

max_forks_repo_name: sophiano/cusvm


text:
#!/usr/bin/env python
import os
import yaml
import numpy as np
import h5py

#ms = range(71, 96+1)
ms = range(98, 126)

os.chdir(os.path.expandvars('$WELTER/sf/'))

for m in ms:

    path_out = 'm{:03d}/'.format(m)
    os.chdir(path_out)

    # Also initialize the nuissance parameter json file.
    os.system('star.py -r 2')
    os.chdir('output/mix_emcee/run02')
    os.system('star.py --initPhi')
    os.chdir(os.path.expandvars('$WELTER/sf/'))

max_stars_repo_name: BrownDwarf/welter

max_issues_repo_name: BrownDwarf/welter

max_forks_repo_name: BrownDwarf/welter


text:
# isort:skip_file
"__st"
# %% [markdown]
"""
# Treex

**Main features**:
* Modules contain their parameters
* Easy transfer learning
* Simple initialization
* No metaclass magic
* No apply method
* No need special versions of `vmap`, `jit`, and friends.

We will showcase each of the above features by creating a very contrived but complete module that will use everything from parameters, states, and random states:
"""

# %%
import typing as tp

import jax.numpy as jnp
import numpy as np

import treex as tx

class NoisyLinear(tx.Module):
    # tree parts are defined by treex annotations
    w: tp.Union[jnp.ndarray, tx.Initializer] = tx.Parameter.node()
    b: tp.Union[jnp.ndarray, tx.Initializer] = tx.Parameter.node()
    key: tp.Union[jnp.ndarray, tx.Initializer] = tx.Rng.node()

    # other annotations are possible but ignored by type
    name: str

    def __init__(self, din, dout):

        # Initializers only expect RNG key
        self.w = tx.Initializer(lambda k: jax.random.uniform(k, shape=(din, dout)))
        self.b = tx.Initializer(lambda k: jax.random.uniform(k, shape=(dout,)))

        # random state is JUST state, we can keep it locally
        self.key = tx.Initializer(lambda k: k)

    def __call__(self, x: np.ndarray) -> np.ndarray:
        assert isinstance(self.key, jnp.ndarray)

        # update state in place
        key, self.key = jax.random.split(self.key, 2)

        # your typical linear operation
        y = jnp.dot(x, self.w) + self.b

        # add noise for fun
        return y + 0.8 * jax.random.normal(key, shape=y.shape)

model = NoisyLinear(1, 1)

print(model)

# %% [markdown]
"""
Initialization is straightforward. The only thing you need to do is to call `init` on your module with a random key:
"""

# %%
import jax

model = model.init(key=42)
print(model)

# %% [markdown]
"""
Now we will be reusing the previous NoisyLinear model, and we will create an optax optimizer that is used to train the model:
"""

# %%
import optax

optimizer = tx.Optimizer(optax.adam(1e-2))

params = model.filter(tx.Parameter)
optimizer = optimizer.init(params)

print(f"{params=}")

# %% [markdown]
"""
Notice that we are already splitting the model into `params` and `states` since we only need to pass the `params` to the optimizer. Next, we will create the loss function, it will take the model parts and the data parts and return the loss plus the new states:
"""
# %%
from functools import partial

@partial(jax.value_and_grad, has_aux=True)
def loss_fn(params: NoisyLinear, model: NoisyLinear, x, y):
    # update params into model
    model = model.merge(params)
    # apply model
    pred_y = model(x)
    # MSE loss
    loss = jnp.mean((y - pred_y) ** 2)
    # return model with state updates
    return loss, model

# %% [markdown]
"""
Notice that we are merging the `params` and `states` into the complete model since we need everything in place to perform the forward pass. Also, we return the updated states from the model. The above steps are required because JAX functional API requires us to be explicit about state management.

**Note**: inside `loss_fn` (wrapped by `value_and_grad`) module can behave like a regular mutable Python object. However, every time it is treated as a pytree a new reference will be created in `jit`, `grad`, `vmap`, etc. It is essential to consider this when using functions like `vmap` inside a module, as JAX will need specific bookkeeping to manage the state correctly.

Next, we will implement the `update` function, it will look indistinguishable from your standard Haiku update, which also separates weights into `params` and `states`:
"""

# %%
@jax.jit
def train_step(model: NoisyLinear, optimizer: tx.Optimizer, x, y):
    # select Parameters
    params = model.filter(tx.Parameter)

    # call loss_fn to get loss, model state, and gradients
    (loss, model), grads = loss_fn(params, model, x, y)

    # apply optax update
    new_params = optimizer.update(grads, params)

    # update new_params
    model = model.merge(new_params)

    return model, optimizer, loss

# %% [markdown]
"""
Before we start training lets get some data:
"""

# %%
import matplotlib.pyplot as plt
import numpy as np

np.random.seed(0)

def get_data(dataset_size: int) -> tp.Tuple[np.ndarray, np.ndarray]:
    x = np.random.normal(size=(dataset_size, 1))
    y = 5 * x - 2 + 0.4 * np.random.normal(size=(dataset_size, 1))
    return x, y

def get_batch(
    data: tp.Tuple[np.ndarray, np.ndarray], batch_size: int
) -> tp.Tuple[np.ndarray, np.ndarray]:
    idx = np.random.choice(len(data[0]), batch_size)
    return jax.tree_map(lambda x: x[idx], data)

data = get_data(1000)

fig = plt.figure()  # __st
plt.scatter(data[0], data[1])
plt.show()
fig  # __st

# %% [markdown]
"""
Finally, we create a simple training loop that performs a few thousand updates and update `params` and `states` back into a single `model` at the end:
"""

# %%
steps = 10_000

for step in range(steps):
    x, y = get_batch(data, batch_size=32)

    model, optimizer, loss = train_step(model, optimizer, x, y)

    if step % 1000 == 0:
        print(f"[{step}] loss = {loss}")

# %% [markdown]
"""
Now lets generate some test data and see how our model performed:
"""

# %%
import matplotlib.pyplot as plt

X_test = np.linspace(data[0].min(), data[0].max(), 100)[:, None]
preds = model(X_test)

fig = plt.figure()  # __st
plt.scatter(data[0], data[1], label="data", color="k")
plt.plot(X_test, preds, label="prediction")
plt.legend()
plt.show()
fig  # __st

# %% [markdown]
"""
We can see that the model has learned the general trend, but because of the `NoisyLinear` modules we have a bit of noise in the predictions.
"""

max_stars_repo_name: ptigwe/treex

max_issues_repo_name: ptigwe/treex

max_forks_repo_name: ptigwe/treex

reason: import numpy,import jax


text:
import datetime
import logging

import matplotlib.pyplot as plt
import numpy as np
import pandas as pd

class AverageDecreaseTermination:

    def __init__(self, N: int, tolerance:  float = 0):
        """ Callback to terminate optimization based the average decrease

        The average decrease over the last N data points is compared to the specified tolerance

        Args:
            N: Number of data points to use
            tolerance: Abort if the average decrease is smaller than the specified tolerance

        """
        self.N = N
        self.tolerance = tolerance
        self.values = []

    def __call__(self, parameters, value) -> bool:
        """
        Returns:
            True if the optimization loop should be aborted
        """
        self.values.append(value)

        if len(self.values) > self.N:
            last_values = self.values[-self.N:]
            pp = np.polyfit(range(self.N), last_values, 1)
            slope = pp[0] / self.N

            logging.debug(f'AverageDecreaseTermination: slope {slope}')
            if slope > self.tolerance:
                return True
        return False

class OptimizerCallback:

    _column_names = ['iteration', 'timestamp', 'residual']

    def __init__(self, show_progress=False, store_data=True):
        """ Class to collected data of optimization procedures """
        self.show_progress = show_progress
        self.store_data = store_data

        self.clear()

    def decorate(self, *args):
        """  Decorate an objective function """
        raise NotImplementedError('')

    @property
    def data(self) -> pd.DataFrame:
        """ Return data gathered by callback """

        df = pd.DataFrame(self._data, columns=self._column_names)

        return df  # self._data

    def _append(self, d):
        self._data.append(d)

    def clear(self):
        self.parameters = []

        self._data = []

    def number_of_evaluations(self) -> int:
        """ Return the number of callback evaluations

        Note: this can differ from the number of objective evaluations
        """
        return len(self.data)

    def optimization_time(self) -> float:
        if len(self.data) > 1:
            delta_t = self.data.iloc[-1]['timestamp']-self.data.iloc[0]['timestamp']
            dt = delta_t.total_seconds()
        else:
            dt = 0
        return dt

    def plot(self, ax=None, **kwargs):
        """ Plot optimization results """
        if ax is None:
            ax = plt.gca()

        self.data.plot('iteration', 'residual', ax=ax, **kwargs)
        dt = self.optimization_time()
        ax.set_title(f'Optimization total time {dt:.2f} [s]')

    def data_callback(self, iteration, parameters, residual):
        """ Callback used to store data """
        if self.store_data:
            self.parameters.append(parameters)

            ts = datetime.datetime.now()  # .isoformat()
            d = (int(iteration), ts, float(residual))
            self._append(d)

    def callback(self, number_evaluations, parameters, value, stepsize, accepted):
        if self.show_progress:
            print(f'#{number_evaluations}, {parameters}, {value}, {stepsize}, {accepted}')

    def qiskit_callback(self, number_evaluations, parameters, value, stepsize, accepted):
        """ Callback method for Qiskit optimizers """
        if self.show_progress:
            print(f'#{number_evaluations}, {parameters}, {value}, {stepsize}, {accepted}')
        self.data_callback(number_evaluations, parameters, value)

    def lmfit_callback(self, parameters, iteration, residual, *args, **kws):
        """ Callback method for lmfit optimizers """
        self.callback(iteration, parameters, residual, 'stepsize', None)
        residual = np.linalg.norm(residual)
        self.data_callback(iteration, parameters, residual)

    def scipy_callback(self, parameters):
        """ Callback method for scipy optimizers """
        number_evaluations = self.number_of_evaluations()
        value = np.NaN
        if self.show_progress:
            print(f'#{number_evaluations}, {parameters}, {value}')
        self.data_callback(number_evaluations, parameters, value)

if __name__ == '__main__':
    import scipy.optimize
    from qiskit.algorithms.optimizers import SPSA

    def rosen(params, a=1, b=100, noise=0):
        """ Rosenbrock function """
        v = (a-params[0])**2+b*(params[1]-params[0]**2)**2
        v += noise*(np.random.rand()-.5)
        return v

    def objective(x):
        return rosen(x, 0.01, 1)

    oc = OptimizerCallback(show_progress=True)

    result = scipy.optimize.minimize(objective, [.5, .9], callback=oc.scipy_callback)

    oc.clear()
    optimizer = SPSA(maxiter=90, callback=oc.qiskit_callback, termination_callback=AverageDecreaseTermination(10))
    point, value, eval = optimizer.optimize(2, objective, initial_point=[0.4, 0.98])
    print(f'point {point}, {value}')
    plt.figure(100)
    plt.clf()
    oc.plot(logy=True)

    print(f'number_of_evaluations {oc.number_of_evaluations()} dt {oc.optimization_time()}')
    from qtt.utilities.tools import measure_time

    with measure_time():
        for ii in range(20):
            oc.data_callback(0, [2, 3], .1)
    with measure_time():
        for ii in range(1000):
            oc.data_callback(0, [2, 3], .1)
    with measure_time():
        for ii in range(20):
            oc.data_callback(0, [2, 3], .1)

    ts = datetime.datetime.now()
    data = [[1, ts, 2], [3, ts, 3]]
    pd.DataFrame(data, columns=['a', 'b', 'c'])

max_stars_repo_name: eendebakpt/qtt

max_issues_repo_name: eendebakpt/qtt

max_forks_repo_name: eendebakpt/qtt

reason: import numpy,import scipy


text:
# coding=utf-8
# Copyright 2019 Google LLC
# 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.

"""Script allowing to play the game by multiple players."""
from __future__ import absolute_import
from __future__ import division
from __future__ import print_function
import os
from baselines.bench import monitor
from baselines import logger
import numpy as np
import pickle
import tensorflow as tf
from tensorflow import Tensor
from tensorflow.keras.layers import Input, Conv2D, ReLU, BatchNormalization, Add, MaxPooling2D, Flatten, Dense
from tensorflow.keras.models import Model
import sonnet as snt
from absl import app
from absl import flags
from absl import logging
from gfootball.env import config
from gfootball.env import football_env
from gfootball.env import wrappers
from gfootball.env import _apply_output_wrappers
from gfootball.env import _process_representation_wrappers
from gfootball.env import observation_preprocessing
from gfootball.examples.models import gfootball_impala_cnn
from gfootball.env.players.ppo2_cnn import Player, ObservationStacker
from gfootball.env.football_action_set import full_action_set
# tf.enable_eager_execution()

# https://www.tensorflow.org/tutorials/quickstart/advanced
# https://towardsdatascience.com/building-a-resnet-in-keras-e8f1322a49ba
# TODO: will unbalanced loading cause problems???

def relu_bn(inputs):
    relu = ReLU()(inputs)
    bn = BatchNormalization()(relu)
    return bn

def residual_block(x, filters, kernel_size=3):
    y = Conv2D(kernel_size=kernel_size, strides=1, filters=filters, padding="same")(x)
    y = relu_bn(y)
    y = Conv2D(kernel_size=kernel_size, strides=1, filters=filters, padding="same")(y)
    out = Add()([x, y])
    out = relu_bn(out)
    return out

def make_gfootball_resnet():

    filters_and_blocks = [(16, 2), (32, 2), (32, 2), (32, 2)]

    inputs = Input(shape=(observation_preprocessing.SMM_HEIGHT, observation_preprocessing.SMM_WIDTH, 16))
    x = inputs

    for n_filters, n_blocks in filters_and_blocks:
        x = Conv2D(kernel_size=3, strides=1, filters=n_filters, padding="same")(x)
        x = MaxPooling2D(pool_size=[3,3], strides=[2,2], padding='same')(x)
        x = relu_bn(x)
        for _ in range(n_blocks):
            x = residual_block(x, n_filters)

    x = ReLU()(x)
    x = Flatten()(x)
    x = Dense(256, activation='relu')(x)
    outputs = Dense(len(full_action_set), activation='softmax')(x)

    model = Model(inputs, outputs)

    return model

def train_on_games(model, start_game, end_game, n_epochs):

    observations = []
    actions = []
    for game_i in range(start_game, end_game):
        with open(f'/Users/stephen/Documents/football/data/observations{game_i}.pkl', 'rb') as f:
            game_obs = pickle.load(f)
            observations.append(game_obs)
        with open(f'/Users/stephen/Documents/football/data/actions{game_i}.pkl', 'rb') as f:
            game_actions = pickle.load(f)
            actions.append(game_actions)

    observations = np.concatenate(observations)
    actions = np.concatenate(actions)

    observations /= 255

    n_total = observations.shape[0]
    n_train = int(0.9 * n_total)
    n_test = n_total - n_train
    train_idxs = np.random.permutation(np.concatenate([np.ones((n_train,), dtype=int),
                                                       np.zeros((n_test,), dtype=int)]))
    train_idxs = train_idxs.astype(bool)

    train_X = observations[train_idxs]
    train_y = actions[train_idxs]
    test_X = observations[1-train_idxs]
    test_y = actions[1-train_idxs]

    loss_weights = 1 / np.bincount(train_y)
    inf_idxs = np.isinf(loss_weights)
    loss_weights[inf_idxs] = 1
    model.compile(optimizer='adam', loss='sparse_categorical_crossentropy', metrics=['accuracy'])
    model.fit(train_X, train_y, batch_size=64, epochs=n_epochs, shuffle=True, validation_data=(test_X, test_y),
              class_weight=loss_weights)

def main(_):

    model = make_gfootball_resnet()
    train_on_games(model, 0, 5, 3)
    train_on_games(model, 5, 10, 3)

    model.save('/Users/stephen/Documents/football/models/behavioral_cloning.h5')

if __name__ == '__main__':
    app.run(main)

max_stars_repo_name: thestephencasper/football

max_issues_repo_name: thestephencasper/football

max_forks_repo_name: thestephencasper/football


text:
from __future__ import print_function

import tensorflow as tf
from tensorflow.contrib import rnn
import numpy as np
import pickle as pickle     # python pkl 文件读写

from tensorflow_LSTMs import MyNetworks

def use():
    train_data = np.array( pickle.load( open('tf_model_lstm/train_seg_data.plk', 'rb') ) )
    train_labels = np.array( pickle.load( open('tf_model_lstm/train_seg_labels.plk', 'rb') ) )

    print( type( train_data ) )
    print( train_data.shape )
    # ndarray

    ''' ********************************************************************************** '''
    batch_size = len(train_data)
    # Training Parameters:
    learning_rate = 0.05  # 0.1 , 0.001

    # Network Parameters
    num_input = 1
    timesteps = 83 # timesteps
    num_hidden = 200 # hidden layer num of features
    num_classes = 4 # MNIST total classes (1-4 digits)

    # 定义神经网络
    train_op, loss_op, prediction, accuracy, X, Y = MyNetworks(learning_rate, num_input, timesteps, num_hidden, num_classes)

    # 声明tf.train.Saver类用于保存/加载模型
    saver = tf.train.Saver()

    # Start training
    # 开始训练
    with tf.Session() as sess:

        # 将固化到硬盘中的Session从保存路径再读取出来
        saver.restore(sess, "tf_model_lstm/model.ckpt")

        batch_x = train_data
        batch_y = train_labels
            # Reshape data to get 28 seq of 28 elements
        batch_x = batch_x.reshape((batch_size, timesteps, num_input))

        print("Testing Accuracy(测试集正确率):", \
            sess.run(accuracy, feed_dict={X: batch_x, Y: batch_y}))

        out = []
        for i in range(batch_size):
            test_data = train_data[i]
            test_label = train_labels[i]
                # Reshape data to get 28 seq of 28 elements
            test_data = test_data.reshape( (1, timesteps, num_input) )

            # print(test_data.shape, test_label)

            # print("分类:", sess.run(prediction, feed_dict={X: test_data}))
            out_tag = sess.run(prediction, feed_dict={X: test_data})
            # 获取矩阵值(概率)最大下标
            j = np.unravel_index( out_tag[0].argmax(), out_tag[0].shape )
            j = j[0]+1
            # print(j)
            out.append(j)

        return out

import matplotlib.pyplot as plt

if __name__ == "__main__":
    print('hello!')
    out = use()

    plt.plot( out )
    plt.show()


text:
import os
from cv2 import cv2
import random
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import h5py
import glob
from PIL import Image

import tensorflow as tf
from keras import backend as K
from keras.models import Model
from keras.layers import Input, Conv2D, MaxPooling2D, Reshape, Bidirectional, LSTM, Dense, Lambda, Activation, BatchNormalization, Dropout
from keras.optimizers import Adam

from keras.models import load_model
from keras.preprocessing import image
import string
from collections import defaultdict
from collections import Counter

from scipy.ndimage import interpolation as inter

from fuzzywuzzy import fuzz
from fuzzywuzzy import process

from Levenshtein import distance
import sys, io
import datetime

def correct_skew(image):
 # Compute rotated bounding box
    thresh = cv2.threshold(image, 0, 255, cv2.THRESH_BINARY + cv2.THRESH_OTSU)[1]
    coords = np.column_stack(np.where(thresh == 0))
    angle = cv2.minAreaRect(coords)[-1]

    if angle < -45:
        angle = -(90 + angle)
    else:
        angle = -angle

    # Rotate image to deskew
    (h, w) = image.shape[:2]
    center = (w // 2, h // 2)
    M = cv2.getRotationMatrix2D(center, angle, 1.0)
    rotated = cv2.warpAffine(image, M, (w, h), flags=cv2.INTER_CUBIC, borderMode=cv2.BORDER_REPLICATE)
    return rotated

def possible_wordsF(fwords, charSet):
    guessed = list()
    scale1 = 0
    scale2= 0
    scale3 = 0
    for word in fwords:
        scale2 = fuzz.ratio(word,charSet)
        if scale2 >= scale1 :
            scale1 = scale2
    for word in fwords:
        scale3 = fuzz.ratio(word,charSet)
        if scale3 == scale1 :
                guessed.append(word)

    print("Content-type:text/plain;charset=utf-8\n\n")
    print("\n Probably : \n")
    print(guessed)

def possible_wordsL(lwords, charSet):
    guessed = list()
    scale1 = len(charSet)
    scale2= 0
    scale3 = 0
    for word in lwords:
        scale2 = distance(word,charSet)
        if len(word)>len(charSet) :
            scale2 = distance(word[0:len(charSet)-1],charSet)
        if scale2 <= scale1 :
               scale1 = scale2
    for word in lwords:
        scale3 = distance(word,charSet)
        if len(word)>len(charSet) :
            scale3 = distance(word[0:len(charSet)-1],charSet)
        if scale3 == scale1 :
                guessed.append(word)
    guessed = '\n'.join(guessed)
    br0 = r'
'
    br1 = r'
'
    br2 = r'
'
    print(br0+" Possible strings : \n"+br2)
    print(br1,guessed,br2)

if __name__ == "__main__":

    # characters = ['0', '1', '2', '3', '4', '5', '6', '7', '8', '9', 'A', 'B', 'C', 'D', 'E', 'F', 'G', 'H', 'I', 'J', 'K',
    #               'L', 'M', 'N', 'O', 'P', 'Q', 'R', 'S', 'T', 'U', 'V', 'W', 'X', 'Y', 'Z', 'a', 'b', 'c', 'd', 'e', 'f',
    #               'g', 'h', 'i', 'j', 'k', 'l', 'm', 'n', 'o', 'p', 'q', 'r', 's', 't', 'u', 'v', 'w', 'x', 'y', 'z']
    characters = [ 'A', 'B', 'C', 'D', 'E', 'F', 'G', 'H', 'I', 'J', 'K',
                    'L', 'M', 'N', 'O', 'P', 'Q', 'R', 'S', 'T', 'U', 'V', 'W', 'X', 'Y', 'Z'
                    ,'a','b','d','e','f','g','h','n','q','r','t']

    # enter input image here
    img = cv2.imread(r'C:\Users\MuktiAS\SKRIPSI\UI\bismillahUI\public_html\obat\obat.png')
    image = cv2.medianBlur(img,7)
    # image = cv2.transpose(image)
    # image = np.asarray(image) # convert from integers to floats
    # image = image.astype('float32') # normalize to the range 0-1 img /= 255.0
    # image = cv2.imread(r'C:\Users\MuktiAS\SKRIPSI\emnist_recog\bw\putra1.jpeg')
    height, width, depth = image.shape

    # resizing the image to find spaces better
    image = cv2.resize(image, dsize=(width * 5, height * 4), interpolation=cv2.INTER_CUBIC)
    # grayscale
    gray = cv2.cvtColor(image, cv2.COLOR_BGR2GRAY)
    # gray = cv2.adaptiveThreshold(gray,255,cv2.ADAPTIVE_THRESH_GAUSSIAN_C,cv2.THRESH_BINARY,11,2)
    # binary
    ret, thresh = cv2.threshold(gray, 127, 255, cv2.THRESH_BINARY_INV)
    #compute a bit-wise inversion so black becomes white and vice versa
    # thresh = np.invert(thresh)

    # dilation
    kernel = np.ones((1, 1), np.uint8)
    img_dilation = cv2.dilate(thresh, kernel, iterations=1)

    # adding GaussianBlur
    gsblur = cv2.GaussianBlur(img_dilation, (5, 5), 0)
    # gsblur = cv2.medianBlur(img_dilation, 9)

    # find contours
    ctrs, hier = cv2.findContours(gsblur.copy(), cv2.RETR_EXTERNAL, cv2.CHAIN_APPROX_SIMPLE)

    m = list()
    # sort contours
    sorted_ctrs = sorted(ctrs, key=lambda ctr: cv2.boundingRect(ctr)[0])
    pchl = list()
    dp = image.copy()
    for i, ctr in enumerate(sorted_ctrs):
        # Get bounding box
        area = cv2.contourArea(ctr)
        # if area > 10000:
        x, y, w, h = cv2.boundingRect(ctr)
        cx = x+w//2
        cy = y+h//2
        cr  = max(w,h)//2
        r = cr+10
        cv2.rectangle(dp, (cx-r, cy-r), (cx+r, cy+r), (90, 0, 255), 9)

    plt.imshow(dp)
    plt.savefig(r'C:\Users\MuktiAS\SKRIPSI\UI\bismillahUI\public_html\obat\obat.png', bbox_inches='tight')

    # Recreate the exact same model, including its weights and the optimizer
    from keras.models import load_model
    from keras.models import model_from_json

    # json_file = open(r'\Users\MuktiAS\SKRIPSI\UI\bismillahUI\app\Http\Controllers\HTR\model95.json')
    # json_file = open(r'\Users\MuktiAS\SKRIPSI\UI\bismillahUI\app\Http\Controllers\HTR\modelDG96In.json')
    # json_file = open(r'\Users\MuktiAS\SKRIPSI\UI\bismillahUI\app\Http\Controllers\HTR\modelDG97InTr.json')
    # json_file = open(r'\Users\MuktiAS\SKRIPSI\UI\bismillahUI\app\Http\Controllers\HTR\modelDG978FlTr.json')
    json_file = open(r'\Users\MuktiAS\SKRIPSI\UI\bismillahUI\app\Http\Controllers\HTR\tesub.json')
    loaded_model_json = json_file.read()
    json_file.close()
    loaded_model = model_from_json(loaded_model_json)
    # loaded_model.load_weights(r'\Users\MuktiAS\SKRIPSI\UI\bismillahUI\app\Http\Controllers\HTR\model95.h5')
    # loaded_model.load_weights(r'\Users\MuktiAS\SKRIPSI\UI\bismillahUI\app\Http\Controllers\HTR\modelDG96In.h5')
    # loaded_model.load_weights(r'\Users\MuktiAS\SKRIPSI\UI\bismillahUI\app\Http\Controllers\HTR\modelDG97InTr.h5')
    # loaded_model.load_weights(r'\Users\MuktiAS\SKRIPSI\UI\bismillahUI\app\Http\Controllers\HTR\modelDG98FlTr.h5')
    loaded_model.load_weights(r'\Users\MuktiAS\SKRIPSI\UI\bismillahUI\app\Http\Controllers\HTR\tesub.h5')
    model = loaded_model
    # model = load_model(r'C:\Users\MuktiAS\SKRIPSI\emnist_recog\dokumentasi\model.h5')
    # model.load_weights(r'C:\Users\MuktiAS\SKRIPSI\emnist_recog\dokumentasi\model_weights.h5')
    # model.summary()

    for i, ctr in enumerate(sorted_ctrs):
        # Get bounding box
        area = cv2.contourArea(ctr)
        # if area > 10000:
        x, y, w, h = cv2.boundingRect(ctr)
        cx = x+w//2
        cy = y+h//2
        cr  = max(w,h)//2
        r = cr+10
        # Getting ROI
        roi = image[(1 if(cy-r) <= 0 else (cy-r)):cy+r, (1 if(cx-r) <= 0 else (cx-r)):cx+r]
        # (1 if(cy-r) <= 0 else (cy-r)):cy + h + 10, (1 if(cx-r) <= 0 else (cx-r)):cx + w + 10
        roi = cv2.resize(roi, dsize=(28, 28), interpolation=cv2.INTER_CUBIC)
        roi = cv2.cvtColor(roi, cv2.COLOR_BGR2GRAY)
        roi1 = cv2.bitwise_not(roi)
        path = 'C:/Users/MuktiAS/SKRIPSI/UI/bismillahUI/public_html/kontur/'
        nama = datetime.datetime.now().strftime("%Y-%m-%d %H-%M-%S-%f")
        # print(nama)
        cv2.imwrite(os.path.join(path , nama+'.png'), roi1)
        # roi = correct_skew(roi)

        roi = np.array(roi)
        t = np.copy(roi)
        t = t / 255.0
        t = 1 - t
        t = t.reshape(1, 28,28)
        m.append(roi)
        pred = np.argmax(model.predict(t.round()), axis=-1)
        pchl.append(pred)

    pcw = list()
    interp = 'bilinear'
    # fig, axs = plt.subplots(nrows=len(sorted_ctrs), sharex=True, figsize=(1, len(sorted_ctrs)))
    for i in range(len(pchl)):
        # print (pchl[i][0])
        pcw.append(characters[pchl[i][0]])
        # axs[i].set_title('-------> predicted letter: ' + characters[pchl[i][0]], x=2.5, y=0.24)
        # axs[i].imshow(m[i], interpolation=interp)

    # plt.show()

    predstring = ''.join(pcw)
    print('Identified String: ' + predstring )
    print('\n')

    tebak = list(predstring)
    # charSet = list(set(tebak))
    # charSet.sort()
    # print(charSet)
    print('\n')

    with open(r'\Users\MuktiAS\SKRIPSI\UI\bismillahUI\app\Http\Controllers\HTR\all_medicine_names.txt') as f:
    # with open(r'\Users\MuktiAS\SKRIPSI\UI\bismillahUI\app\Http\Controllers\HTR\300_obat.txt') as f:
        input = [line.strip() for line in f]
    # print(input)

    # possible_wordsF(input, tebak)
    possible_wordsL(input, predstring)

reason: import numpy,from scipy