text:
#!/usr/bin/python
import matplotlib.pyplot as plt
import os
import sys
sys.path.insert(0, '/Users/paytonrodman/athena/vis/python')
import csv
import numpy as np
import scipy.stats
from ast import literal_eval
import argparse

def main(**kwargs):
    # directory containing data
    problem  = args.prob_id
    root_dir = '/Users/paytonrodman/athena-sim/'
    data_dir = root_dir + problem + '/'
    os.chdir(data_dir)
    time = []
    Bav = []
    Bcc1_0sh = []
    Bcc1_1sh = []
    Bcc1_2sh = []
    Bcc1_3sh = []
    with open('B_strength_with_time.csv', newline='\n') as csv_file:
        csv_reader = csv.reader(csv_file, delimiter='\t')
        next(csv_reader, None) # skip header
        for row in csv_reader:
            t = float(row[0])
            bav = row[1]
            bcc1_0sh = row[2]
            bcc1_1sh = row[3]
            bcc1_2sh = row[4]
            bcc1_3sh = row[5]

            time.append(t)
            Bav.append(literal_eval(bav))
            Bcc1_0sh.append(float(bcc1_0sh))
            Bcc1_1sh.append(literal_eval(bcc1_1sh))
            Bcc1_2sh.append(literal_eval(bcc1_2sh))
            Bcc1_3sh.append(literal_eval(bcc1_3sh))

    time, Bav, Bcc1_0sh, Bcc1_1sh, Bcc1_2sh, Bcc1_3sh = zip(*sorted(zip(time, Bav, Bcc1_0sh, Bcc1_1sh, Bcc1_2sh, Bcc1_3sh)))

    fig, axs = plt.subplots(2)
    Bav = np.abs(np.array(Bav))
    axs[1].plot(time,Bav[:,0],label=r'$B_{\Phi}$')
    axs[1].plot(time,Bcc1_0sh,label=r'$B_{\rm midplane}$')
    axs[1].plot(time,np.array(Bcc1_1sh)[:,0],label=r'$B_{1H/r}$')
    axs[1].plot(time,np.array(Bcc1_2sh)[:,0],label=r'$B_{2H/r}$')
    axs[1].plot(time,np.array(Bcc1_3sh)[:,0],label=r'$B_{3H/r}$')
    axs[1].set_xlabel(r'time ($GM/c^3$)')
    axs[1].set_ylabel(r'$\langle B\rangle_{r<100r_g}$ in disk (lower)')
    axs[1].legend()

    axs[0].plot(time,Bav[:,1],label=r'$B_{\Phi}$')
    axs[0].plot(time,Bcc1_0sh,label=r'$B_{\rm midplane}$')
    axs[0].plot(time,np.array(Bcc1_1sh)[:,1],label=r'$B_{1H/r}$')
    axs[0].plot(time,np.array(Bcc1_2sh)[:,1],label=r'$B_{2H/r}$')
    axs[0].plot(time,np.array(Bcc1_3sh)[:,1],label=r'$B_{3H/r}$')
    axs[0].set_xlabel(r'time ($GM/c^3$)')
    axs[0].set_ylabel(r'$\langle B\rangle_{r<100r_g}$ in disk (upper)')
    axs[0].legend()

    for ax in fig.get_axes():
        ax.label_outer()
        ax.grid(b=True, which='major', color='#666666', linestyle='-', alpha=0.5)
        ax.minorticks_on()
        ax.ticklabel_format(axis="x", style="sci", scilimits=(0,0))
        ax.grid(b=True, which='minor', color='#999999', linestyle='-', alpha=0.2)

    plt.tight_layout()
    plt.savefig(data_dir + 'B_strength' + '.png', dpi=1200)
    plt.close()

    if args.logy:
        fig = plt.figure()
        ax = fig.add_subplot(111)
        ax.semilogy(time,beta)
        ax.grid()
        ax.set_xlabel(r'time ($GM/c^3$)')
        ax.set_ylabel(r'average $\beta$ in disk')
        plt.ticklabel_format(axis="x", style="sci", scilimits=(0,0))
        plt.grid(b=True, which='major', color='#666666', linestyle='-', alpha=0.5)
        plt.minorticks_on()
        plt.grid(b=True, which='minor', color='#999999', linestyle='-', alpha=0.2)
        plt.savefig(data_dir + 'logy_beta' + '.png', dpi=1200)
        plt.close()

def mean_confidence_interval(data, confidence=0.95):
    a = 1.0 * np.array(data)
    n = len(a)
    m, se = np.nanmean(a), scipy.stats.sem(a, nan_policy='omit')
    h = se * scipy.stats.t.ppf((1 + confidence) / 2., n-1)
    return m-h, m+h

# Execute main function
if __name__ == '__main__':
    parser = argparse.ArgumentParser(description='Plot B value over time.')
    parser.add_argument('prob_id',
                        help='base name of the data being analysed, e.g. inflow_var or disk_base')
    parser.add_argument('--logy',
                        action='store_true',
                        help='plot logy version')
    args = parser.parse_args()

    main(**vars(args))

max_stars_repo_name: paytonrodman/athena-analysis

max_issues_repo_name: paytonrodman/athena-analysis

max_forks_repo_name: paytonrodman/athena-analysis

reason: import numpy,import scipy


text:
#!/usr/bin/python3
# -*- coding: utf-8 -*-
'''
@File: demo.py
@Author:kong
@Time: 2020年01月21日09时40分
@Description:
'''
import os
import sys
module_path = os.path.abspath(os.path.join('..'))
if module_path not in sys.path:
    sys.path.append(module_path)
from matplotlib import pyplot as plt
from data import VOCDetection, VOC_ROOT, VOCAnnotationTransform
import torch
import torch.nn as nn
import torch.backends.cudnn as cudnn
from torch.autograd import Variable
import numpy as np
import cv2
if torch.cuda.is_available():
    torch.set_default_tensor_type('torch.cuda.FloatTensor')
from data import VOC_CLASSES as labels
from ssd import build_ssd

image_path = './test/example.jpg'
weight_path = './weights/ssd300_VOC_100000.pth'
model_input = 300

net = build_ssd('test', model_input, 21)    # initialize SSD
net.load_weights(weight_path)
image = cv2.imread(image_path, cv2.IMREAD_COLOR)  # uncomment if dataset not downloaded
rgb_image = cv2.cvtColor(image, cv2.COLOR_BGR2RGB)
x = cv2.resize(image, (model_input, model_input)).astype(np.float32)
x -= (104.0, 117.0, 123.0)
x = x.astype(np.float32)
x = x[:, :, ::-1].copy()
x = torch.from_numpy(x).permute(2, 0, 1)

xx = Variable(x.unsqueeze(0))     # wrap tensor in Variable
if torch.cuda.is_available():
    xx = xx.cuda()
y = net(xx)

top_k=10
detections = y.data
# scale each detection back up to the image
scale = torch.Tensor(rgb_image.shape[1::-1]).repeat(2)    #4个尺度的缩放系数
for i in range(detections.size(1)):          #遍历num_class
    j = 0
    while detections[0,i,j,0] >= 0.2:
        score = detections[0,i,j,0]
        label_name = labels[i-1]
        display_txt = '%s: %.2f'%(label_name, score)
        pt = (detections[0,i,j,1:]*scale).cpu().numpy()
        j+=1
        image = cv2.rectangle(image,(pt[0],pt[1]),(pt[2],pt[3]),(255,0,0),2)
        image = cv2.putText(image,display_txt,(pt[2],pt[1]),cv2.FONT_HERSHEY_COMPLEX,1,(255,0,0),2)
cv2.imwrite('./test/resut.jpg',image)

ext: py

reason: import numpy


text:
import numpy as np
import tensorflow as tf
from tensorflow.contrib.framework.python.ops import add_arg_scope, arg_scope
import math

from matplotlib import pyplot as plt

import tensorflow.contrib.distributions as tfd

def model_arg_scope(**kwargs):
    """Create new counter and apply arg scope to all arg scoped nn
    operations."""
    counters = {}
    return arg_scope(
        [conv2d, deconv2d, residual_block, dense, activate], counters=counters, **kwargs
    )

def make_model(name, template, **kwargs):
    """Create model with fixed kwargs."""
    run = lambda *args, **kw: template(
        *args, **dict((k, v) for kws in (kw, kwargs) for k, v in kws.items())
    )
    if tf.executing_eagerly():
        return tf.make_template(name, run)
    return tf.make_template(name, run, unique_name=name)

def int_shape(x):
    """ short for x.shape.as_list() """
    return x.shape.as_list()

def get_name(layer_name, counters):
    """ utlity for keeping track of layer names """
    if not layer_name in counters:
        counters[layer_name] = 0
    name = layer_name + "_" + str(counters[layer_name])
    counters[layer_name] += 1
    return name

def apply_partwise(input_, func):
    """
    Applies function func on all parts separately.
    Parts are in channel 3.
    The input is reshaped to map the parts to the batch axis and then the function is applied
    Parameters
    ----------
    input_ : tensor
        [b, h, w, parts, features]
    func : callable
        a NN function to apply to each part individually
    Returns
    -------
        [b, out_h, out_w, parts, out_features]
    """

    b, h, w, parts, f = input_.shape.as_list()

    # transpose [b, h, w, part, features] --> [part, b, h, w, features]
    perm = [3, 0, 1, 2, 4]
    x = tf.transpose(input_, perm=perm)
    # reshape [part, b, h, w, features] --> [part * b, h, w, features]
    x = tf.reshape(x, (b * parts, h, w, f))

    y = func(x)

    _, h_out, w_out, c_out = y.shape.as_list()
    # reshape [part * b, h_out, w_out, c_out] --> [part, b, h_out, w_out, c_out]
    out = tf.reshape(y, (parts, b, h_out, w_out, c_out))
    # transpose back [part, b, h_out, w_out, c_out] --> [b, h_out, w_out, part, c_out]
    inv_perm = [1, 2, 3, 0, 4]
    out = tf.transpose(out, perm=inv_perm)
    return out

@add_arg_scope
def partwise_conv2d(
    x,
    num_filters,
    filter_size=[3, 3],
    stride=[1, 1],
    pad="SAME",
    init_scale=1.0,
    counters={},
    init=False,
    initdist="uniform",
    **kwargs
):
    """
    input: [b, h, w, parts, features]
    Each part (channel 3) has is its own bias and scale
    Uses 3D convolution internally to prevent tf.transpose

    Examples
    --------

    import tensorflow as tf
    tf.enable_eager_execution()

    from pylab import *
    from skimage import data
    import numpy as np
    import math
    im = data.astronaut()
    im = im.astype(np.float32) / 255
    H, W, D = im.shape

    b = 1
    parts = 5
    out_features = 1
    features = tf.reshape(im, (b, H, W, 1, D))
    features = tf.concat([features] * parts, axis=3)

    out = partwise_conv2d(features, out_features, init=False, part_wise=True, initdist="debug")

    fig, ax = plt.subplots(1, 1, figsize=(20, 20))
    a = np.hstack([np.squeeze(out[..., p, :]) for p in range(parts)])
    ax.imshow(a, cmap=plt.cm.gray) # this should render the astronaut in 5 different shades of gray

    """
    num_filters = int(num_filters)
    name = get_name("conv2d", counters)
    with tf.variable_scope(name):
        in_channels = x.shape.as_list()[4]
        in_parts = x.shape.as_list()[3]
        fan_in = in_channels * filter_size[0] * filter_size[1]
        stdv = math.sqrt(1.0 / fan_in)
        part_stdv = math.sqrt(1.0 / in_parts)
        if initdist == "uniform":
            V_initializer = tf.random_uniform_initializer(minval=-stdv, maxval=stdv)
            b_initializer = tf.random_uniform_initializer(minval=-stdv, maxval=stdv)
            v_part_initializer = tf.random_uniform_initializer(
                minval=-part_stdv, maxval=part_stdv
            )
            b_part_initializer = tf.random_uniform_initializer(
                minval=-part_stdv, maxval=part_stdv
            )
        elif initdist == "normal":
            V_initializer = tf.random_normal_initializer(stddev=stdv)
            b_initializer = tf.random_normal_initializer(stddev=stdv)
            v_part_initializer = tf.random_normal_initializer(stddev=part_stdv)
            b_part_initializer = tf.random_normal_initializer(stddev=part_stdv)
        elif initdist == "debug":
            pass
        else:
            raise ValueError(initdist)
        if not initdist == "debug":
            V = tf.get_variable(
                "V",
                filter_size + [1, in_channels, num_filters],
                initializer=V_initializer,
                dtype=tf.float32,
            )
            b = tf.get_variable(
                "b",
                [1, 1, 1, 1, num_filters],
                initializer=b_initializer,
                dtype=tf.float32,
            )
            V_part = tf.get_variable(
                "V_part",
                [1, 1, 1, in_parts, 1],
                initializer=v_part_initializer,
                dtype=tf.float32,
            )
            b_part = tf.get_variable(
                "b_part",
                [1, 1, 1, in_parts, 1],
                initializer=b_part_initializer,
                dtype=tf.float32,
            )
        else:
            V = (
                tf.reshape([1.0, 2.0, 1.0], (3, 1, 1, 1, 1))
                / 4
                * tf.reshape([1.0, 2.0, 1.0], (1, 3, 1, 1, 1))
                / 4
                * tf.reshape([1.0, 1.0, 0.0], (1, 1, 1, 3, 1))
                / 2
            )
            b = tf.zeros([1, 1, 1, 1, num_filters], dtype=tf.float32)
            V_part = tf.reshape(
                tf.cast(tf.linspace(0.0, 1.0, in_parts), dtype=tf.float32),
                [1, 1, 1, in_parts, 1],
            )
            b_part = tf.zeros([1, 1, 1, in_parts, 1], dtype=tf.float32)
        x = tf.nn.conv3d(x, V, strides=[1] + stride + [1] + [1], padding="SAME")
        x *= V_part
        x += b + b_part
        return x

def _conv2d(
    x,
    num_filters,
    filter_size=[3, 3],
    stride=[1, 1],
    pad="SAME",
    init_scale=1.0,
    counters={},
    init=False,
    **kwargs
):
    num_filters = int(num_filters)
    strides = [1] + stride + [1]
    name = get_name("conv2d", counters)
    initdist = "uniform"
    with tf.variable_scope(name):
        in_channels = int(x.get_shape()[-1])
        fan_in = in_channels * filter_size[0] * filter_size[1]
        stdv = math.sqrt(1.0 / fan_in)
        if initdist == "uniform":
            V_initializer = tf.random_uniform_initializer(minval=-stdv, maxval=stdv)
            b_initializer = tf.random_uniform_initializer(minval=-stdv, maxval=stdv)
        elif initdist == "normal":
            V_initializer = tf.random_normal_initializer(stddev=stdv)
            b_initializer = tf.random_normal_initializer(stddev=stdv)
        else:
            raise ValueError(initdist)
        V = tf.get_variable(
            "V",
            filter_size + [in_channels, num_filters],
            initializer=V_initializer,
            dtype=tf.float32,
        )
        b = tf.get_variable(
            "b", [num_filters], initializer=b_initializer, dtype=tf.float32
        )
        if init:
            tmp = tf.nn.conv2d(x, V, [1] + stride + [1], pad) + tf.reshape(
                b, [1, 1, 1, num_filters]
            )
            mean, var = tf.nn.moments(tmp, [0, 1, 2])
            scaler = 1.0 / tf.sqrt(var + 1e-6)
            V = tf.assign(V, V * scaler)
            b = tf.assign(b, -mean * scaler)
        x = tf.nn.conv2d(x, V, [1] + stride + [1], pad) + tf.reshape(
            b, [1, 1, 1, num_filters]
        )
        return x

@add_arg_scope
def conv2d(
    x,
    num_filters,
    filter_size=[3, 3],
    stride=[1, 1],
    pad="SAME",
    init_scale=1.0,
    counters={},
    init=False,
    part_wise=False,
    coords=False,
    **kwargs
):
    """
    A 2D convolution.

    Parameters
    ----------
    x: tensor
        input tensor [N, H, W, C]
    num_filters: int
        number of feature maps
    filter_size: list of `ints`
        filter size in x, y
    stride: list of `ints`
        stride in x, y
    pad: A `string` from: `"SAME", "VALID"`.
      The type of padding algorithm to use.
        option 'SAME'
    init_scale
    counters
    init
    part_wise: bool
        if True, the input has to be [N, H, W, Parts, C]. The convolution will get an additional scale and bias per part
    coords: bool
        if True, will use coordConv (2018ACS_liuIntriguingFailingConvolutionalNeuralNetworks)
    **kwargs

    Returns
    -------
    tensor
        convolved input

    """
    if coords:
        x = add_coordinates(x)
    if part_wise:
        out = partwise_conv2d(
            x,
            num_filters,
            filter_size=filter_size,
            stride=stride,
            pad=pad,
            init_scale=init_scale,
            counters=counters,
            init=init,
            **kwargs
        )
    else:
        out = _conv2d(
            x,
            num_filters,
            filter_size=filter_size,
            stride=stride,
            pad=pad,
            init_scale=init_scale,
            counters=counters,
            init=init,
            **kwargs
        )

    return out

@add_arg_scope
def dense(x, num_units, init_scale=1.0, counters={}, init=False, **kwargs):
    """ fully connected layer """
    name = get_name("dense", counters)
    initdist = "uniform"
    with tf.variable_scope(name):
        in_channels = int(x.get_shape()[-1])
        fan_in = in_channels
        stdv = math.sqrt(1.0 / fan_in)
        if initdist == "uniform":
            V_initializer = tf.random_uniform_initializer(minval=-stdv, maxval=stdv)
            b_initializer = tf.random_uniform_initializer(minval=-stdv, maxval=stdv)
        elif initdist == "normal":
            V_initializer = tf.random_normal_initializer(stddev=stdv)
            b_initializer = tf.random_normal_initializer(stddev=stdv)
        else:
            raise ValueError(initdist)
        V = tf.get_variable(
            "V", [in_channels, num_units], initializer=V_initializer, dtype=tf.float32
        )
        b = tf.get_variable(
            "b", [num_units], initializer=b_initializer, dtype=tf.float32
        )
        if init:
            tmp = tf.matmul(x, V) + tf.reshape(b, [1, num_units])
            mean, var = tf.nn.moments(tmp, [0])
            scaler = 1.0 / tf.sqrt(var + 1e-6)
            V = tf.assign(V, V * scaler)
            b = tf.assign(b, -mean * scaler)
        x = tf.matmul(x, V) + tf.reshape(b, [1, num_units])
        return x

@add_arg_scope
def activate(x, activation, **kwargs):
    """
    Activation unit
    Parameters
    ----------
    x: tensor
        input tensor
    activation:
        A `string` from: `"elu", "relu", "leaky_relu", "softplus"`.
    kwargs

    Returns
    -------

    """
    if activation == None:
        return x
    elif activation == "elu":
        return tf.nn.elu(x)
    elif activation == "relu":
        return tf.nn.relu(x)
    elif activation == "leaky_relu":
        return tf.nn.leaky_relu(x)
    elif activation == "softplus":
        return tf.nn.softplus(x)
    else:
        raise NotImplemented(activation)

def nin(x, num_units):
    """ a network in network layer (1x1 CONV) """
    return conv2d(x, num_units, filter_size=[1, 1])

def downsample(x, num_units):
    """
    Downsampling by stride 2 convolution

    equivalent to
    x = conv2d(x, num_units, stride = [2, 2])

    Parameters
    ----------
    x: tensor
        input
    num_units:
        number of feature map in the output

    Returns
    -------

    """
    return conv2d(x, num_units, stride=[2, 2])

def upsample(x, num_units, method="subpixel"):
    """
    2D upsampling layer.

    Parameters
    ----------
    x: tensor
        input
    num_units:
        number of feature maps in the output
    method:
        upsampling method. A `string` from: `"conv_transposed", "nearest_neighbor", "linear", "subpixel"`
        Subpixel means that every upsampled pixel gets its own filter.

    Returns
    -------
        upsampled input

    """
    xs = x.shape.as_list()
    if method == "conv_transposed":
        return deconv2d(x, num_units, stride=[2, 2])
    elif method == "subpixel":
        x = conv2d(x, 4 * num_units)
        x = tf.depth_to_space(x, 2)
        return x
    elif method == "nearest_neighbor":
        bs, h, w, c = x.shape.as_list()
        x = tf.image.resize_images(
            x, [2 * h, 2 * w], tf.image.ResizeMethod.NEAREST_NEIGHBOR
        )
        return x
    elif method == "linear":
        bs, h, w, c = xs[:4]
        x = tf.image.resize_images(x, [2 * h, 2 * w], tf.image.ResizeMethod.BILINEAR)
        return x
    else:
        raise NotImplemented(method)

@add_arg_scope
def partwise_deconv2d(
    x,
    num_filters,
    filter_size=[3, 3],
    stride=[1, 1],
    pad="SAME",
    init_scale=1.0,
    counters={},
    init=False,
    **kwargs
):
    """ transposed convolutional layer """
    num_filters = int(num_filters)
    name = get_name("deconv2d", counters)
    xs = x.shape.as_list()
    strides = [1] + stride + [1]
    in_parts = xs[3]
    part_stdv = math.sqrt(1.0 / in_parts)
    v_part_initializer = tf.random_uniform_initializer(
        minval=-part_stdv, maxval=part_stdv
    )
    b_part_initializer = tf.random_uniform_initializer(
        minval=-part_stdv, maxval=part_stdv
    )
    if pad == "SAME":
        target_shape = [
            xs[0] * in_parts,
            xs[1] * stride[0],
            xs[2] * stride[1],
            num_filters,
        ]
    else:
        target_shape = [
            xs[0] * in_parts,
            xs[1] * stride[0] + filter_size[0] - 1,
            xs[2] * stride[1] + filter_size[1] - 1,
            num_filters,
        ]
    with tf.variable_scope(name):
        V = tf.get_variable(
            "V",
            filter_size + [num_filters, xs[-1]],
            tf.float32,
            tf.random_normal_initializer(0, 0.05),
        )
        g = tf.get_variable(
            "g",
            [num_filters],
            dtype=tf.float32,
            initializer=tf.constant_initializer(1.0),
        )
        b = tf.get_variable(
            "b",
            [num_filters],
            dtype=tf.float32,
            initializer=tf.constant_initializer(0.0),
        )
        V_part = tf.get_variable(
            "V_part", [in_parts], initializer=v_part_initializer, dtype=tf.float32
        )
        b_part = tf.get_variable(
            "b_part", [in_parts], initializer=b_part_initializer, dtype=tf.float32
        )
        V_norm = tf.nn.l2_normalize(V, [0, 1, 3])

        def part_conv_func(x_):
            x_ = tf.nn.conv2d_transpose(
                x_, V_norm, target_shape, [1] + stride + [1], pad
            )
            x_ = tf.reshape(g, [1, 1, 1, num_filters]) * x_ + tf.reshape(
                b, [1, 1, 1, num_filters]
            )
            return x_

        if init:
            mean, var = tf.nn.moments(x, [0, 1, 2])
            g = tf.assign(g, init_scale / tf.sqrt(var + 1e-10))
            b = tf.assign(b, -mean * g)
        x = apply_partwise(x, part_conv_func)
        x = x * tf.reshape(V_part, (1, 1, 1, in_parts, 1)) + tf.reshape(
            b_part, (1, 1, 1, in_parts, 1)
        )
        return x

@add_arg_scope
def _deconv2d(
    x,
    num_filters,
    filter_size=[3, 3],
    stride=[1, 1],
    pad="SAME",
    init_scale=1.0,
    counters={},
    init=False,
    **kwargs
):
    """ transposed convolutional layer """
    num_filters = int(num_filters)
    name = get_name("deconv2d", counters)
    xs = x.shape.as_list()
    strides = [1] + stride + [1]
    if pad == "SAME":
        target_shape = [xs[0], xs[1] * stride[0], xs[2] * stride[1], num_filters]
    else:
        target_shape = [
            xs[0],
            xs[1] * stride[0] + filter_size[0] - 1,
            xs[2] * stride[1] + filter_size[1] - 1,
            num_filters,
        ]
    with tf.variable_scope(name):
        V = tf.get_variable(
            "V",
            filter_size + [num_filters, xs[-1]],
            tf.float32,
            tf.random_normal_initializer(0, 0.05),
        )
        g = tf.get_variable(
            "g",
            [num_filters],
            dtype=tf.float32,
            initializer=tf.constant_initializer(1.0),
        )
        b = tf.get_variable(
            "b",
            [num_filters],
            dtype=tf.float32,
            initializer=tf.constant_initializer(0.0),
        )

        V_norm = tf.nn.l2_normalize(V, [0, 1, 3])
        x = tf.nn.conv2d_transpose(x, V_norm, target_shape, [1] + stride + [1], pad)
        if init:
            mean, var = tf.nn.moments(x, [0, 1, 2])
            g = tf.assign(g, init_scale / tf.sqrt(var + 1e-10))
            b = tf.assign(b, -mean * g)
        x = tf.reshape(g, [1, 1, 1, num_filters]) * x + tf.reshape(
            b, [1, 1, 1, num_filters]
        )
        return x

@add_arg_scope
def deconv2d(
    x,
    num_filters,
    filter_size=[3, 3],
    stride=[1, 1],
    pad="SAME",
    init_scale=1.0,
    counters={},
    init=False,
    part_wise=False,
    coords=False,
    **kwargs
):
    """
        coords: if True, will use coordConv (2018ACS_liuIntriguingFailingConvolutionalNeuralNetworks)
    """
    if coords:
        x = add_coordinates(x)
    if part_wise:
        out = partwise_deconv2d(
            x,
            num_filters,
            filter_size=filter_size,
            stride=stride,
            pad=pad,
            init_scale=init_scale,
            counters=counters,
            init=init,
            **kwargs
        )
    else:
        out = _deconv2d(
            x,
            num_filters,
            filter_size=filter_size,
            stride=stride,
            pad=pad,
            init_scale=init_scale,
            counters=counters,
            init=init,
            **kwargs
        )
    return out

@add_arg_scope
def residual_block(x, skipin=None, conv=conv2d, init=False, dropout_p=0.0, **kwargs):
    """
    slight variation of original residual block.

    Parameters
    ----------
    x: tensor
        Incoming tensor
    skipin: tensor
        Incomming tensor from skip connection, if any
    conv: callable
        which convolution function to use for the resodual
    init
    dropout_p: float
        dropout probability, if any
    kwargs

    Returns
    -------
    output : x + residual

    Examples

    """
    xs = int_shape(x)
    num_filters = xs[-1]

    residual = x
    if skipin is not None:
        skipin = nin(activate(skipin), num_filters)
        residual = tf.concat([residual, skipin], axis=-1)
    residual = activate(residual)
    residual = tf.nn.dropout(residual, keep_prob=1.0 - dropout_p)
    residual = conv(residual, num_filters)

    return x + residual

def flatten(x):
    """
    returns a flat version of x --> [N, -1]
    Parameters
    ----------
    x: tensor

    Returns
    -------

    """
    _shape = x.shape.as_list()
    return tf.reshape(x, (_shape[0], -1))

def mask2rgb(mask):
    """
    Convert tensor with masks [N, H, W, C] to an RGB tensor [N, H, W, 3]
    using argmax over channels.
    Parameters
    ----------
    mask: ndarray
        an array of shape [N, H, W, C]
    Returns:
        RGB visualization in shape [N, H, W, 3]
    -------

    """
    n_parts = mask.shape.as_list()[3]
    maxmask = tf.argmax(mask, axis=3)
    hotmask = tf.one_hot(maxmask, depth=n_parts)
    hotmask = tf.expand_dims(hotmask, 4)
    colors = make_mask_colors(n_parts)
    colors = (colors - 0.5) * 2
    colors = tf.to_float(colors)
    colors = tf.expand_dims(colors, 0)
    colors = tf.expand_dims(colors, 0)
    colors = tf.expand_dims(colors, 0)
    vis_mask = hotmask * colors
    vis_mask = tf.reduce_sum(vis_mask, axis=3)
    return vis_mask

def np_one_hot(targets, n_classes):
    """
    numpy equivalent of tf.one_hot
    returns targets as one hot matrix

    Parameters
    ----------
    targets: ndarray
        array of target classes
    n_classes: int
        how many classes there are overall
    Returns: ndarray
        one-hot array with shape [n, n_classes]
    -------

    """
    res = np.eye(n_classes)[np.array(targets).reshape(-1)]
    return res.reshape(list(targets.shape) + [n_classes])

def np_to_float(x):
    """ cast x to float32 """
    return x.astype(np.float32)

def np_mask2rgb(mask):
    """
    numpy equivalent of @mask2rgb
    convert tensor with masks [N, H, W, C] to an RGB tensor [N, H, W, 3]
    using argmax over channels.
    Parameters
    ----------
    mask: ndarray
        an array of shape [N, H, W, C]
    Returns:
        RGB visualization in shape [N, H, W, 3]
    -------

    """
    n_parts = mask.shape[3]
    maxmask = np.argmax(mask, axis=3)
    hotmask = np_one_hot(maxmask, n_classes=n_parts)
    hotmask = np.expand_dims(hotmask, 4)
    colors = make_mask_colors(n_parts)
    colors = (colors - 0.5) * 2
    colors = np_to_float(colors)
    colors = np.expand_dims(colors, 0)
    colors = np.expand_dims(colors, 0)
    colors = np.expand_dims(colors, 0)
    vis_mask = hotmask * colors
    vis_mask = np.sum(vis_mask, axis=3)
    return vis_mask

def make_mask_colors(n_parts, cmap=plt.cm.inferno):
    """
    make a color array using the specified colormap for n_parts classes
    Parameters
    ----------
    n_parts: int
        how many classes there are in the mask
    cmap:
        matplotlib colormap handle

    Returns
    -------
    colors: ndarray
        an array with shape [n_parts, 3] representing colors in the range [0, 1].

    """
    colors = cmap(np.linspace(0, 1, n_parts), alpha=False, bytes=False)[:, :3]
    return colors

def hourglass_model(
    x,
    config,
    extra_resnets,
    n_out=3,
    activation="relu",
    upsample_method="subpixel",
    coords=False,
):
    """
    A U-net or hourglass style image-to-image model with skip-connections

    Parameters
    ----------
    x : tensor
        input tensor to unet
    config : list
        a list of ints specifying the number of feature maps on each scale of the unet in the downsampling path
        for the upsampling path, the list will be reversed
        For example [32, 64] will use 32 channels on scale 0 (without downsampling) and 64 channels on scale 1
        once downsampled).
    extra_resnets : int
        how many extra res blocks to use at the bottleneck
    n_out : int
        number of final output feature maps of the unet. 3 for RGB
    activation : str
        a string specifying the activation function to use. See @activate for options.
    upsample_method : list of str or str
        a str specifying the upsampling method or a list of str specifying the upsampling method for each scale individually.
        See @upsample for possible options.
    coords : True
        if coord conv should be used.

    Returns
    -------

    Examples
    --------

    tf.enable_eager_execution()
    x = tf.ones((1, 128, 128, 3))
    config = [32, 64]
    extra_resnets = 0
    upsample_method = "subpixel"
    activation = "leaky_relu"
    coords = False

    unet = make_model("unet", hourglass_model, config=config, extra_resnets= extra_resnets, upsample_method=upsample_method, activation=activation)
    y = unet(x)

    # plotting the output should look random because we did not train anything
    im = np.concatenate([x, y], axis=1)
    plt.imshow(np.squeeze(im))

    """
    with model_arg_scope(activation=activation, coords=coords):
        hs = list()

        h = conv2d(x, config[0])
        h = residual_block(h)

        for nf in config[1:]:
            h = downsample(h, nf)
            h = residual_block(h)
            hs.append(h)

        for _ in range(extra_resnets):
            h = residual_block(h)

        for i, nf in enumerate(config[-2::-1]):
            h = residual_block(h, skipin=hs[-(i + 1)])
            h = upsample(h, nf, method=upsample_method)

        h = residual_block(h)
        h = conv2d(h, n_out)
        return h

def make_ema(init_value, value, decay=0.99):
    """
    apply exponential moving average to variable

    Parameters
    ----------
    init_value: float
        initial value for moving average variable
    value: variable
        tf variable to apply update ops on
    decay: float
        decay parameter

    Returns
    -------
    avg_value : variable with exponential moving average
    update_ema: tensorflow update operation for exponential moving average

    Examples
    --------

    # usage within edflow Trainer.make_loss_ops. Apply EMA to discriminator accuracy
    avg_acc, update_ema = make_ema(0.5, dis_accuracy, decay)
    self.update_ops.append(update_ema)
    self.log_ops["dis_acc"] = avg_acc

    """
    decay = tf.constant(decay, dtype=tf.float32)
    avg_value = tf.Variable(init_value, dtype=tf.float32, trainable=False)
    update_ema = tf.assign(
        avg_value, decay * avg_value + (1.0 - decay) * tf.cast(value, tf.float32)
    )
    return avg_value, update_ema

def add_coordinates(input_tensor, with_r=False):
    """
    Given an input_tensor, adds 2 channelw ith x and y coordinates to the feature maps.
    This was introduced in coordConv (2018ACS_liuIntriguingFailingConvolutionalNeuralNetworks).
    Parameters
    ----------
    input_tensor: tensor
        Tensor of shape [N, H, W, C]
    with_r : bool
        if True, euclidian radius will also be added as channel

    Returns
    ret : input_tensor concatenated with x and y coordinates and maybe euclidian distance.
    -------

    """
    assert len(input_tensor.shape.as_list()) == 4
    bs = tf.shape(input_tensor)[0]
    x_dim, y_dim = input_tensor.shape.as_list()[1:3]
    xx_ones = tf.ones([bs, x_dim], dtype=tf.int32)
    xx_ones = tf.expand_dims(xx_ones, -1)
    xx_range = tf.tile(tf.expand_dims(tf.range(y_dim), 0), [bs, 1])
    xx_range = tf.expand_dims(xx_range, 1)

    xx_channel = tf.matmul(xx_ones, xx_range)
    xx_channel = tf.expand_dims(xx_channel, -1)

    yy_ones = tf.ones([bs, y_dim], dtype=tf.int32)
    yy_ones = tf.expand_dims(yy_ones, 1)
    yy_range = tf.tile(tf.expand_dims(tf.range(x_dim), 0), [bs, 1])
    yy_range = tf.expand_dims(yy_range, -1)

    yy_channel = tf.matmul(yy_range, yy_ones)
    yy_channel = tf.expand_dims(yy_channel, -1)

    xx_channel = tf.cast(xx_channel, "float32") / max(1, x_dim - 1)
    yy_channel = tf.cast(yy_channel, "float32") / max(1, y_dim - 1)
    xx_channel = xx_channel * 2 - 1
    yy_channel = yy_channel * 2 - 1

    ret = tf.concat([input_tensor, xx_channel, yy_channel], axis=-1)
    if with_r:
        rr = tf.sqrt(tf.square(xx_channel) + tf.square(yy_channel))
        ret = tf.concat([ret, rr], axis=-1)
    return ret

def probs_to_mu_L(
    probs, scaling_factor, inv=True
):  # todo maybe exponential map induces to much certainty ! low values basically ignored and only high values count!
    """Calculate mean and covariance (cholesky decomposition of covariance) for each channel of probs
    tensor of keypoint probabilites [bn, h, w, n_kp]
    mean calculated on a grid of scale [-1, 1]

    Parameters
    ----------
    probs: tensor
        tensor of shape [b, h, w, k] where each channel along axis 3 is interpreted as an unnormalized probability density.
    scaling_factor : tensor
        tensor of shape [b, 1, 1, k] representing normalizing the normalizing constant of the density
    inv: bool
        if True, returns covariance matrix of density. Else returns inverse of covariance matrix aka precision matrix

    Returns
    -------
    mu : tensor
        tensor of shape [b, k, 2] representing partwise mean coordinates of x and y for each item in the batch
    L : tensor
        tensor of shape [b, k, 2, 2] representing partwise cholesky decomposition of covariance
         matrix for each item in the batch.

    Example
    -------

    .. code-block:: python

        from matplotlib import pyplot as plt
        tf.enable_eager_execution()
        import numpy as np
        import tensorflow.contrib.distributions as tfd

        _means = [-0.5, 0, 0.5]
        means = tf.ones((3, 1, 2), dtype=tf.float32) * np.array(_means).reshape((3, 1, 1))
        means = tf.concat([means, means, means[::-1, ...]], axis=1)
        means = tf.reshape(means, (-1, 2))

        var_ = 0.1
        rho = 0.5
        cov = [[var_, rho * var_], [rho * var_, var_]]
        scale = tf.cholesky(cov)
        scale = tf.stack([scale] * 3, axis=0)
        scale = tf.stack([scale] * 3, axis=0)
        scale = tf.reshape(scale, (-1, 2, 2))

        mvn = tfd.MultivariateNormalTriL(
            loc=means,
            scale_tril=scale)

        h = 100
        w = 100
        y_t = tf.tile(tf.reshape(tf.linspace(-1., 1., h), [h, 1]), [1, w])
        x_t = tf.tile(tf.reshape(tf.linspace(-1., 1., w), [1, w]), [h, 1])
        y_t = tf.expand_dims(y_t, axis=-1)
        x_t = tf.expand_dims(x_t, axis=-1)
        meshgrid = tf.concat([y_t, x_t], axis=-1)
        meshgrid = tf.expand_dims(meshgrid, 0)
        meshgrid = tf.expand_dims(meshgrid, 3)  # 1, h, w, 1, 2

        blob = mvn.prob(meshgrid)
        blob = tf.reshape(blob, (100, 100, 3, 3))
        blob = tf.transpose(blob, perm=[2, 0, 1, 3])

        norm_const = np.sum(blob, axis=(1, 2), keepdims=True)
        mu, L = nn.probs_to_mu_L(blob / norm_const, 1, inv=False)

        bn, h, w, nk = blob.get_shape().as_list()
        estimated_blob = nn.tf_hm(h, w, mu, L)

        fig, ax = plt.subplots(2, 3, figsize=(9, 6))
        for b in range(len(_means)):
            ax[0, b].imshow(np.squeeze(blob[b, ...]))
            ax[0, b].set_title("target_blobs")
            ax[0, b].set_axis_off()

        for b in range(len(_means)):
            ax[1, b].imshow(np.squeeze(estimated_blob[b, ...]))
            ax[1, b].set_title("estimated_blobs")
            ax[1, b].set_axis_off()

    """
    (
        bn,
        h,
        w,
        nk,
    ) = (
        probs.get_shape().as_list()
    )  # todo instead of calulating sequrity measure from amplitude one could alternativly calculate it by letting the network predict a extra paremeter also one could do
    y_t = tf.tile(tf.reshape(tf.linspace(-1.0, 1.0, h), [h, 1]), [1, w])
    x_t = tf.tile(tf.reshape(tf.linspace(-1.0, 1.0, w), [1, w]), [h, 1])
    y_t = tf.expand_dims(y_t, axis=-1)
    x_t = tf.expand_dims(x_t, axis=-1)
    meshgrid = tf.concat([y_t, x_t], axis=-1)

    mu = tf.einsum("ijl,aijk->akl", meshgrid, probs)
    mu_out_prod = tf.einsum(
        "akm,akn->akmn", mu, mu
    )  # todo incosisntent ordereing of mu! compare with cross_V2

    mesh_out_prod = tf.einsum(
        "ijm,ijn->ijmn", meshgrid, meshgrid
    )  # todo efficient (expand_dims)
    stddev = tf.einsum("ijmn,aijk->akmn", mesh_out_prod, probs) - mu_out_prod

    a_sq = stddev[:, :, 0, 0]
    a_b = stddev[:, :, 0, 1]
    b_sq_add_c_sq = stddev[:, :, 1, 1]
    eps = 1e-12  # todo clean magic

    a = tf.sqrt(
        a_sq + eps
    )  # Σ = L L^T Prec = Σ^-1  = L^T^-1 * L^-1  ->looking for L^-1 but first L = [[a, 0], [b, c]
    b = a_b / (a + eps)
    c = tf.sqrt(b_sq_add_c_sq - b ** 2 + eps)
    z = tf.zeros_like(a)

    if inv:
        det = tf.expand_dims(tf.expand_dims(a * c, axis=-1), axis=-1)
        row_1 = tf.expand_dims(
            tf.concat(
                [tf.expand_dims(c, axis=-1), tf.expand_dims(z, axis=-1)], axis=-1
            ),
            axis=-2,
        )
        row_2 = tf.expand_dims(
            tf.concat(
                [tf.expand_dims(-b, axis=-1), tf.expand_dims(a, axis=-1)], axis=-1
            ),
            axis=-2,
        )

        L_inv = (
            scaling_factor / (det + eps) * tf.concat([row_1, row_2], axis=-2)
        )  # L^⁻1 = 1/(ac)* [[c, 0], [-b, a]
        return mu, L_inv
    else:
        row_1 = tf.expand_dims(
            tf.concat(
                [tf.expand_dims(a, axis=-1), tf.expand_dims(z, axis=-1)], axis=-1
            ),
            axis=-2,
        )
        row_2 = tf.expand_dims(
            tf.concat(
                [tf.expand_dims(b, axis=-1), tf.expand_dims(c, axis=-1)], axis=-1
            ),
            axis=-2,
        )

        L = scaling_factor * tf.concat([row_1, row_2], axis=-2)  # just L
        return mu, L

def probs_to_mu_sigma(probs):
    """Calculate mean and covariance matrix for each channel of spatial probability maps
    Mean and covariance are caluclated on a grid of scale [-1, 1]

    Parameters
    ----------
    probs: tensor
        tensor of shape [N, H, W, C] where each channel along axis 3 is interpreted as a probability density.

    Returns
    -------
    mu : tensor
        tensor of shape [N, C, 2] representing partwise mean coordinates of x and y for each item in the batch
    sigma : tensor
        tensor of shape [N, C, 2, 2] representing covariance matrix matrix for each item in the batch

    Example
    -------

        mu, sigma = nn.probs_to_mu_sigma(spatial_probability_maps)
    """
    bn, h, w, nk = probs.get_shape().as_list()
    y_t = tf.tile(tf.reshape(tf.linspace(-1.0, 1.0, h), [h, 1]), [1, w])
    x_t = tf.tile(tf.reshape(tf.linspace(-1.0, 1.0, w), [1, w]), [h, 1])
    y_t = tf.expand_dims(y_t, axis=-1)
    x_t = tf.expand_dims(x_t, axis=-1)
    meshgrid = tf.concat([y_t, x_t], axis=-1)

    mu = tf.einsum("ijl,aijk->akl", meshgrid, probs)
    mu_out_prod = tf.einsum("akm,akn->akmn", mu, mu)

    mesh_out_prod = tf.einsum("ijm,ijn->ijmn", meshgrid, meshgrid)
    sigma = tf.einsum("ijmn,aijk->akmn", mesh_out_prod, probs) - mu_out_prod
    return mu, sigma

def tf_hm(h, w, mu, L):
    """
    Returns Gaussian densitiy function based on μ and L for each batch index and part
    L is the cholesky decomposition of the covariance matrix : Σ = L L^T

    Parameters
    ----------
    h : int
        heigh ot output map
    w : int
        width of output map
    mu : tensor
        mean of gaussian part and batch item. Shape [b, p, 2]. Mean in range [-1, 1] with respect to height and width
    L : tensor
        cholesky decomposition of covariance matrix for each batch item and part. Shape [b, p, 2, 2]
    order:

    Returns
    -------
    density : tensor
        gaussian blob for each part and batch idx. Shape [b, h, w, p]

    Example
    -------

    .. code-block:: python

        from matplotlib import pyplot as plt
        tf.enable_eager_execution()
        import numpy as np
        import tensorflow as tf
        import tensorflow.contrib.distributions as tfd

        # create Target Blobs
        _means = [-0.5, 0, 0.5]
        means = tf.ones((3, 1, 2), dtype=tf.float32) * np.array(_means).reshape((3, 1, 1))
        means = tf.concat([means, means, means[::-1, ...]], axis=1)
        means = tf.reshape(means, (-1, 2))

        var_ = 0.1
        rho = 0.5
        cov = [[var_, rho * var_],
               [rho * var_, var_]]
        scale = tf.cholesky(cov)
        scale = tf.stack([scale] * 3, axis=0)
        scale = tf.stack([scale] * 3, axis=0)
        scale = tf.reshape(scale, (-1, 2, 2))

        mvn = tfd.MultivariateNormalTriL(
            loc=means,
            scale_tril=scale)

        h = 100
        w = 100
        y_t = tf.tile(tf.reshape(tf.linspace(-1., 1., h), [h, 1]), [1, w])
        x_t = tf.tile(tf.reshape(tf.linspace(-1., 1., w), [1, w]), [h, 1])
        y_t = tf.expand_dims(y_t, axis=-1)
        x_t = tf.expand_dims(x_t, axis=-1)
        meshgrid = tf.concat([y_t, x_t], axis=-1)
        meshgrid = tf.expand_dims(meshgrid, 0)
        meshgrid = tf.expand_dims(meshgrid, 3)  # 1, h, w, 1, 2

        blob = mvn.prob(meshgrid)
        blob = tf.reshape(blob, (100, 100, 3, 3))
        blob = tf.transpose(blob, perm=[2, 0, 1, 3])

        # Estimate mean and L
        norm_const = np.sum(blob, axis=(1, 2), keepdims=True)
        mu, L = nn.probs_to_mu_L(blob / norm_const, 1, inv=False)

        bn, h, w, nk = blob.get_shape().as_list()

        # Estimate blob based on mu and L
        estimated_blob = nn.tf_hm(h, w, mu, L)

        # plot
        fig, ax = plt.subplots(2, 3, figsize=(9, 6))
        for b in range(len(_means)):
            ax[0, b].imshow(np.squeeze(blob[b, ...]))
            ax[0, b].set_title("target_blobs")
            ax[0, b].set_axis_off()

        for b in range(len(_means)):
            ax[1, b].imshow(np.squeeze(estimated_blob[b, ...]))
            ax[1, b].set_title("estimated_blobs")
            ax[1, b].set_axis_off()

    """

    assert len(mu.get_shape().as_list()) == 3
    assert len(L.get_shape().as_list()) == 4
    assert mu.get_shape().as_list()[-1] == 2
    assert L.get_shape().as_list()[-1] == 2
    assert L.get_shape().as_list()[-2] == 2

    b, p, _ = mu.get_shape().as_list()
    mu = tf.reshape(mu, (b * p, 2))
    L = tf.reshape(L, (b * p, 2, 2))

    mvn = tfd.MultivariateNormalTriL(loc=mu, scale_tril=L)
    y_t = tf.tile(tf.reshape(tf.linspace(-1.0, 1.0, h), [h, 1]), [1, w])
    x_t = tf.tile(tf.reshape(tf.linspace(-1.0, 1.0, w), [1, w]), [h, 1])
    y_t = tf.expand_dims(y_t, axis=-1)
    x_t = tf.expand_dims(x_t, axis=-1)
    meshgrid = tf.concat([y_t, x_t], axis=-1)
    meshgrid = tf.expand_dims(meshgrid, 0)
    meshgrid = tf.expand_dims(meshgrid, 3)  # 1, h, w, 1, 2

    probs = mvn.prob(meshgrid)
    probs = tf.reshape(probs, (h, w, b, p))
    probs = tf.transpose(probs, perm=[2, 0, 1, 3])  # move part axis to the back
    return probs

if __name__ == "__main__":
    tf.enable_eager_execution()
    x = tf.ones((1, 128, 128, 3))
    config = [32, 64]
    extra_resnets = 0
    upsample_method = "subpixel"
    activation = "leaky_relu"
    coords = False

    unet = make_model(
        "unet",
        hourglass_model,
        config=config,
        extra_resnets=extra_resnets,
        upsample_method=upsample_method,
        activation=activation,
    )
    y = unet(x)

    im = np.concatenate([x, y], axis=1)
    plt.imshow(np.squeeze(im))

max_stars_repo_name: theRealSuperMario/edflow

max_issues_repo_name: theRealSuperMario/edflow

max_forks_repo_name: ArWeHei/edflow


text:
import warnings
from typing import List

import numpy as np
import tensorflow as tf
from tfgraph.utils.math.vector_norm import VectorNorm

from tfgraph.algorithms.pagerank import Transition
from tfgraph.graph import Graph
from tfgraph.utils.callbacks.update_edge_listener import UpdateEdgeListener
from tfgraph.utils.math.convergence_criterion import ConvergenceCriterion
from tfgraph.utils.tensorflow_object import TensorFlowObject
from tfgraph.utils.utils import Utils

class PageRank(TensorFlowObject, UpdateEdgeListener):
  """ PageRank base class.

  This class model the PageRank algorithm as Abstract Class containing all
  methods that the heir classes need to implements. Also, this class provides
  a set of attributes that helps to implement the algorithm.

  The PageRank algorithm calculates the rank of each vertex in a graph based on
  the relational structure from them and giving more importance to the vertices
  that connects with edges to vertices with very high in-degree recursively.

  This class depends on the TensorFlow library, so it's necessary to install it
  to properly work.

  Attributes:
    sess (:obj:`tf.Session`): This attribute represents the session that runs
      the TensorFlow operations.
    name (str): This attribute represents the name of the object in TensorFlow's
      op Graph.
    G (:obj:`tfgraph.Graph`): The graph on witch it will be calculated the
      algorithm. It will be treated as Directed Weighted Graph.
    beta (float): The reset probability of the random walks, i.e. the
      probability that a user that surfs the graph an decides to jump to another
      vertex not connected to the current.
    T (:obj:`tfgraph.Transition`): The transition matrix that provides the
      probability distribution relative to the walk to another node of the graph.
    v (:obj:`tf.Variable`): The stationary distribution vector. It contains the
      normalized probability to stay in each vertex of the graph. So represents
      the PageRank ranking of the graph.
    writer (:obj:`tf.summary.FileWriter`): This attribute represents a
      TensorFlow's Writer, that is used to obtain stats.
    is_sparse (bool): Use sparse Tensors if it's set to True. Not
      implemented yet.

  """

  def __init__(self, sess: tf.Session, name: str, beta: float,
               T: Transition, writer: tf.summary.FileWriter = None,
               is_sparse: bool = False) -> None:
    """ The constructor of the class.

    This method initializes all the attributes needed to compute the PageRank
    of the graph.

    Args:
      sess (:obj:`tf.Session`): This attribute represents the session that runs
        the TensorFlow operations.
      name (str): This attribute represents the name of the object in
        TensorFlow's op Graph.
      beta (float): The reset probability of the random walks, i.e. the
        probability that a user that surfs the graph an decides to jump to
        another vertex not connected to the current.
      T (:obj:`tfgraph.Transition`): The transition matrix that provides the
        probability distribution relative to the walk to another node of the
        graph.
      v (:obj:`tf.Variable`): The stationary distribution vector. It contains
        the normalized probability to stay in each vertex of the graph. So
        represents the PageRank ranking of the graph.
      writer (:obj:`tf.summary.FileWriter`): This attribute represents a
        TensorFlow's Writer, that is used to obtain stats.
      is_sparse (bool): Use sparse Tensors if it's set to True. Not
        implemented yet.

    """
    TensorFlowObject.__init__(self, sess, name, writer=writer,
                              is_sparse=is_sparse)
    UpdateEdgeListener.__init__(self)

    self.beta = beta
    self.T = T
    self.T.attach(self)
    self.v = tf.Variable(tf.fill([1, self.T.G.n], tf.pow(self.T.G.n_tf, -1)),
                         name=self.T.G.name + "_" + self.name + "_v")
    self.run_tf(tf.variables_initializer([self.v]))

  def error_vector_compare_tf(self, other_pr: 'PageRank',
                              k: int = -1) -> tf.Tensor:
    """ The comparison method between two PageRank algorithm results.

    This method compares the `self` PageRank with another one passed as
    parameter of the function. The comparison is based on the difference of the
    Norm One of each `v` vector.

    The method also provides a `k` parameter as option to base the comparison
    only the `k` better ranked vertices.

    Args:
      other_pr (:obj:`tfgraph.PageRank`): Another PageRank object to compare the
        resulting ranking.
      k (int, optional): An additional parameter that allows to base the
        comparison only on the `k` better vertices. Not implemented yet.

    Returns:
      (:obj:`tf.Tensor`): A `tf.Tensor` with 0-D shape, that represents the
        difference between the two rankings using the Norm One.

    Todo:
      * Implement ranking based only on the `k` better ranked vertices.

    """
    if 0 < k < self.T.G.n - 1:
      if 0 < k < self.T.G.n - 1:
        warnings.warn('k-best error comparison not implemented yet')

    return tf.reshape(
      VectorNorm.ONE(tf.subtract(self.v, other_pr.v)), [])

  def error_vector_compare_np(self, other_pr: 'PageRank',
                              k: int = -1) -> np.ndarray:
    """ The comparison method between two PageRank algorithm results.

    This method compares the `self` PageRank with another one passed as
    parameter of the function. The comparison is based on the difference of the
    Norm One of each `v` vector.

    The method also provides a `k` parameter as option to base the comparison
    only the `k` better ranked vertices.

    Args:
      other_pr (:obj:`tfgraph.PageRank`): Another PageRank object to compare the
        resulting ranking.
      k (int, optional): An additional parameter that allows to base the
        comparison only on the `k` better vertices. Not implemented yet.

    Returns:
      (:obj:`np.ndarray`): A `np.ndarray` with 0-D shape, that represents the
        difference between the two rankings using the Norm One.

    """
    return self.run_tf(self.error_vector_compare_tf(other_pr, k))

  def pagerank_vector_tf(self, convergence: float = 1.0, steps: int = 0,
                         topics: List[int] = None,
                         topics_decrement: bool = False,
                         c_criterion=ConvergenceCriterion.ONE) -> tf.Tensor:
    """ The Method that runs the PageRank algorithm

    This method generates a TensorFlow graph of operations needed to calculate
    the PageRank Algorithm and sets to it different parameters passed as
    parameters.

    This method acts as interface between the algorithm and the external
    classes, so it contains a set of parameters that in some implementations of
    PageRank algorithms will not be needed. All the parameters is defined as
    optional for this reason.

    Args:
      convergence (float, optional): A float between 0 and 1 that represents
        the convergence rate that allowed to finish the iterative
        implementations of the algorithm to accept the solution. It has more
        preference than the `steps` parameter. Default to `1.0`.
      steps (int, optional): A positive integer that sets the number of
        iterations that the iterative implementations will run the algorithm
        until finish. It has less preference than the `convergence` parameter.
        Default to `0`.
      topics (:obj:`list` of :obj:`int`, optional): A list of integers that
        represent the set of vertex where the random jumps arrives. If this
        parameter is used, the uniform distribution over all vertices of the
        random jumps will be modified to jump only to this vertex set. Default
        to `None`.
      topics_decrement (bool, optional): If topics is not None and
        topics_decrement is `True` the topics will be casted to 0-Index. Default `
        to False`.
      c_criterion (:obj:`function`, optional): The function used to calculate if
        the Convergence Criterion of the iterative implementations is reached.
        Default to `tfgraph.ConvergenceCriterion.ONE`.

    Returns:
      (:obj:`tf.Tensor`): A 1-D `tf.Tensor` of [n] shape, where `n` is the
        cardinality of the graph vertex set. It contains the normalized rank of
        vertex `i` at position `i`.

    """
    if topics_decrement is True and topics is not None:
      topics = [item - 1 for item in topics]

    if 0.0 < convergence < 1.0:
      return self._pr_convergence_tf(convergence, topics=topics,
                                     c_criterion=c_criterion)
    elif steps > 0:
      return self._pr_steps_tf(steps, topics=topics)
    else:
      return self._pr_exact_tf(topics=topics)

  def pagerank_vector_np(self, convergence: float = 1.0, steps: int = 0,
                         topics: List[int] = None,
                         c_criterion=ConvergenceCriterion.ONE) -> np.ndarray:
    """ The Method that runs the PageRank algorithm

    This method returns a Numpy Array that contains the result of running the
    PageRank algorithm customized by the parameters passed to it.

    This method acts as interface between the algorithm and the external
    classes, so it contains a set of parameters that in some implementations of
    PageRank algorithms will not be needed. All the parameters is defined as
    optional for this reason.

    Args:
      convergence (float, optional): A float between 0 and 1 that represents
        the convergence rate that allowed to finish the iterative
        implementations of the algorithm to accept the solution. It has more
        preference than the `steps` parameter. Default to `1.0`.
      steps (int, optional): A positive integer that sets the number of
        iterations that the iterative implementations will run the algorithm
        until finish. It has less preference than the `convergence` parameter.
        Default to `0`.
      topics (:obj:`list` of :obj:`int`, optional): A list of integers that
        represent the set of vertex where the random jumps arrives. If this
        parameter is used, the uniform distribution over all vertices of the
        random jumps will be modified to jump only to this vertex set. Default
        to `None`.
      c_criterion (:obj:`function`, optional): The function used to calculate if
        the Convergence Criterion of the iterative implementations is reached.
        Default to `tfgraph.ConvergenceCriterion.ONE`.

    Returns:
      (:obj:`np.ndarray`): A 1-D `np.ndarray` of [n] shape, where `n` is the
        cardinality of the graph vertex set. It contains the normalized rank of
        vertex `i` at position `i`.

    """
    return self.run_tf(
      self.pagerank_vector_tf(convergence, steps, topics,
                              c_criterion))

  def ranks_np(self, convergence: float = 1.0, steps: int = 0,
               topics: List[int] = None,
               topics_decrement: bool = False) -> np.ndarray:
    """ Generates a ranked version of PageRank results.

    This method returns the PageRank ranking of the graph sorted by the position
    of each vertex in the rank. So it generates a 2-D matrix with shape [n,2]
    where n is the cardinality of the vertex set of the graph, and at the first
    column it contains the index of vertex and the second column contains it
    normalized rank. The `i` row is referred to the vertex with `i` position in
    the rank.

    Args:
      convergence (float, optional): A float between 0 and 1 that represents
        the convergence rate that allowed to finish the iterative
        implementations of the algorithm to accept the solution. It has more
        preference than the `steps` parameter. Default to `1.0`.
      steps (int, optional): A positive integer that sets the number of
        iterations that the iterative implementations will run the algorithm
        until finish. It has less preference than the `convergence` parameter.
        Default to `0`.
      topics (:obj:`list` of :obj:`int`, optional): A list of integers that
        represent the set of vertex where the random jumps arrives. If this
        parameter is used, the uniform distribution over all vertices of the
        random jumps will be modified to jump only to this vertex set. Default
        to `None`.
      topics_decrement (bool, optional): If topics is not None and
        topics_decrement is `True` the topics will be casted to 0-Index. Default `
        to False`.

    Returns:
      (:obj:`np.ndarray`): A 2-D `np.ndarray` than represents a sorted PageRank
        ranking of the graph.

    """
    self.pagerank_vector_tf(convergence, steps, topics, topics_decrement)
    ranks = tf.map_fn(
      lambda x: [x, tf.gather(tf.reshape(self.v, [self.T.G.n]), x)],
      tf.transpose(
        tf.py_func(Utils.ranked, [tf.scalar_mul(-1, self.v)], tf.int64)),
      dtype=[tf.int64, tf.float32])
    return np.concatenate(self.run_tf(ranks), axis=1)

  def _pr_convergence_tf(self, convergence: float, topics: List[int] = None,
                         c_criterion=ConvergenceCriterion.ONE) -> tf.Tensor:
    """ Abstract method to implement a iterative version of PageRank until
      convergence rate.

    This method runs the PageRank algorithm in iterative fashion a undetermined
    number of times bounded by the `convergence` rate and the 'c_criterion'
    criterion.

    Args:
      convergence (float): A float between 0 and 1 that represents
        the convergence rate that allowed to finish the iterative
        implementations of the algorithm to accept the solution. Default to
        `1.0`.
      topics (:obj:`list` of :obj:`int`, optional): A list of integers that
        represent the set of vertex where the random jumps arrives. If this
        parameter is used, the uniform distribution over all vertices of the
        random jumps will be modified to jump only to this vertex set. Default
        to `None`.
      c_criterion (:obj:`function`, optional): The function used to calculate if
        the Convergence Criterion of the iterative implementations is reached.
        Default to `tfgraph.ConvergenceCriterion.ONE`.

    Returns:
      (:obj:`tf.Tensor`): A 1-D `tf.Tensor` of [n] shape, where `n` is the
        cardinality of the graph vertex set. It contains the normalized rank of
        vertex `i` at position `i`.

    """
    raise NotImplementedError(
      'subclasses must override page_rank_until_convergence()!')

  def _pr_steps_tf(self, steps: int, topics: List[int] = None) -> tf.Tensor:
    """ Abstract method to implement a iterative version of PageRank with fixed
      steps.

    This method runs the PageRank algorithm in iterative fashion a fixed number
    of times bounded by the `steps` parameter.

    Args:
      steps (int): A positive integer that sets the number of
        iterations that the iterative implementations will run the algorithm
        until finish. Default to `0`.
      topics (:obj:`list` of :obj:`int`, optional): A list of integers that
        represent the set of vertex where the random jumps arrives. If this
        parameter is used, the uniform distribution over all vertices of the
        random jumps will be modified to jump only to this vertex set. Default
        to `None`.

    Returns:
      (:obj:`tf.Tensor`): A 1-D `tf.Tensor` of [n] shape, where `n` is the
        cardinality of the graph vertex set. It contains the normalized rank of
        vertex `i` at position `i`.

    """
    raise NotImplementedError(
      'subclasses must override page_rank_until_steps()!')

  def _pr_exact_tf(self, topics: List[int] = None) -> tf.Tensor:
    """ Abstract method to implement a exact version of PageRank.

    This method calculates the PageRank of the graph in exact mode.

    Args:
      topics (:obj:`list` of :obj:`int`, optional): A list of integers that
        represent the set of vertex where the random jumps arrives. If this
        parameter is used, the uniform distribution over all vertices of the
        random jumps will be modified to jump only to this vertex set. Default
        to `None`.

    Returns:
      (:obj:`tf.Tensor`): A 1-D `tf.Tensor` of [n] shape, where `n` is the
        cardinality of the graph vertex set. It contains the normalized rank of
        vertex `i` at position `i`.

    """
    raise NotImplementedError(
      'subclasses must override page_rank_exact()!')

  def update_edge(self, edge: np.ndarray, change: float) -> None:
    """ The callback to receive notifications about edge changes in the graph.

    This method is called from the Graph when an addition or deletion is
    produced on the edge set. So probably is necessary to recompute the PageRank
    ranking.

    Args:
      edge (:obj:`np.ndarray`): A 1-D `np.ndarray` that represents the edge that
        changes in the graph, where `edge[0]` is the source vertex, and
        `edge[1]` the destination vertex.
      change (float): The variation of the edge weight. If the final value is
        0.0 then the edge is removed.

    Returns:
      This method returns nothing.

    """
    raise NotImplementedError(
      'subclasses must override update_edge()!')

max_stars_repo_name: tfgraph/tfgraph

max_issues_repo_name: tfgraph/tfgraph

max_forks_repo_name: tfgraph/tfgraph


text:
#!/usr/bin/env python3

from io import StringIO
from itertools import groupby, permutations, product
from math import cos, sin, radians
from sys import argv

import numpy as np

INTERSECTION_SIZE = 12
MATCHING_DIST_COUNT = (INTERSECTION_SIZE * (INTERSECTION_SIZE - 1)) / 2

SIN_90 = sin(radians(90))
COS_90 = cos(radians(90))
X_ROT = np.matrix([[1, 0, 0], [0, COS_90, -SIN_90], [0, SIN_90, COS_90]])
Y_ROT = np.matrix([[COS_90, 0, SIN_90], [0, 1, 0], [-SIN_90, 0, COS_90]])
Z_ROT = np.matrix([[COS_90, -SIN_90, 0], [SIN_90, COS_90, 0], [0, 0, 1]])

def point_diff(a, b):
    return tuple(coord_a - coord_b for (coord_a, coord_b) in zip(a, b))

def sq_dist(a, b):
    return sum(coord**2 for coord in point_diff(a, b))

def manhattan_dist(a, b):
    return sum(abs(coord) for coord in point_diff(a, b))

def group_len(grouper):
    return sum(1 for _ in grouper)

def get_dist_freqs(points):
    sq_distances = [sq_dist(a, b) for a, b in permutations(points, 2)]
    return {d: group_len(g) for d, g in groupby(sorted(sq_distances))}

def get_diffs(ref_points, points):
    pos_diffs = [point_diff(a, b) for a, b in product(ref_points, points)]
    return [(d, group_len(g)) for d, g in groupby(sorted(pos_diffs))]

def rotate_point_cloud(points, rot):
    rotated = np.rint(points.copy() @ rot.T)
    return np.array(rotated, dtype=np.int32)

def align_point_cloud(ref_points, points):
    for _ in range(4):
        for _ in range(4):
            for _ in range(4):
                diff = max(get_diffs(ref_points, points), key=lambda it: it[1])
                if diff[1] >= INTERSECTION_SIZE:
                    return diff[0], points

                points = rotate_point_cloud(points, Y_ROT)
            points = rotate_point_cloud(points, X_ROT)
        points = rotate_point_cloud(points, Z_ROT)

    raise Exception("Alignment not found")

class PointCloud:
    def __init__(self, points):
        self.points = points
        self.dist_freqs = get_dist_freqs(points)
        self.scanner_pos = (0, 0, 0)

    def does_merge_with(self, other_cloud):
        self_freq, other_freq = self.dist_freqs, other_cloud.dist_freqs
        common_dists = self_freq.keys() & other_freq.keys()
        if not common_dists:
            return False

        common_dist_count = sum(min(self_freq[k], other_freq[k]) for k in common_dists)
        return common_dist_count >= MATCHING_DIST_COUNT

    def merge_into(self, scanner):
        diff, points = align_point_cloud(scanner.points, self.points)
        self.scanner_pos = diff

        translated = np.add(points.copy(), np.array(diff))
        all_points = np.unique(np.append(scanner.points, translated, axis=0), axis=0)

        scanner.points = all_points
        scanner.dist_freqs = get_dist_freqs(all_points)

def parse_scanner_data(section):
    section_buffer = StringIO(section)
    points = np.loadtxt(section_buffer, skiprows=1, delimiter=",", dtype=np.int32)
    return PointCloud(points)

def merge_point_clouds(scanners):
    target, *rest = scanners
    done = []
    while rest:
        cloud = rest.pop(0)
        if cloud.does_merge_with(target):
            cloud.merge_into(target)
            done.append(cloud)
        else:
            rest.append(cloud)

    return target, done

def find_max_manhattan_dist(scanners):
    return max(
        manhattan_dist(a.scanner_pos, b.scanner_pos)
        for a, b in permutations(scanners, 2)
    )

with open(argv[1]) as input_file:
    sections = input_file.read().split("\n\n")
    point_clouds = [parse_scanner_data(section) for section in sections]

scanner_0, _ = merge_point_clouds(point_clouds)
print(f"Part 1: {scanner_0.points.shape[0]}")

max_distance = find_max_manhattan_dist(point_clouds)
print(f"Part 2: {max_distance}")


text:
import warnings
from collections import OrderedDict

import numpy
from scipy.optimize import root

from dolo.numeric.misc import MyJacobian

def find_deterministic_equilibrium(model, constraints=None,
                                   return_jacobian=False):
    '''
    Finds the steady state calibration.

    Taking the value of parameters as given, finds the values for endogenous
    variables consistent with the deterministic steady-state.

    This function requires the specification of the first order equations.

    Parameters
    ----------
    model: NumericModel
        an `(f,g)` compliant model
    constraints: dict
        a dictionaries with forced values.
        Use it to set shocks to non-zero values or to add additional
        constraints in order to avoid unit roots.

    Returns:
    --------
    OrderedDict:
        calibration dictionary (i.e. endogenous variables and parameters by
        type)
    '''

    f = model.functions['arbitrage']
    g = model.functions['transition']

    s0 = model.calibration['states']
    x0 = model.calibration['controls']
    p = model.calibration['parameters']

    if 'shocks' in model.calibration:
        e0 = model.calibration['shocks'].copy()
    else:
        e0 = numpy.zeros(len(model.symbols['shocks']))

    n_e = len(e0)

    z = numpy.concatenate([s0, x0, e0])

    symbs = model.symbols['states'] + model.symbols['controls']
    addcons_ind = []
    addcons_val = []

    if constraints is None:
        constraints = dict()

    for k in constraints:
        if k in symbs:
            i = symbs.index(k)
            addcons_ind.append(i)
            addcons_val.append(constraints[k])
        elif k in model.symbols['shocks']:
            i = model.symbols['shocks'].index(k)
            e0[i] = constraints[k]
        else:
            raise Exception(
                "Invalid symbol '{}' for steady_state constraint".format(k))

    def fobj(z):
        s = z[:len(s0)]
        x = z[len(s0):-n_e]
        e = z[-n_e:]

        S = g(s, x, e, p)
        r = f(s, x, e, s, x, p)
        d_e = e - e0
        d_sx = z[addcons_ind] - addcons_val
        res = numpy.concatenate([S - s, r, d_e, d_sx])
        return res

    jac = MyJacobian(fobj)(z)

    if return_jacobian:
        return jac

    rank = numpy.linalg.matrix_rank(jac)

    if rank < len(z):
        msg = """\
        There are {} equilibrium variables to find, but the jacobian \
        matrix is only of rank {}. The solution is indeterminate."""
        warnings.warn(msg.format(len(z), rank))

    sol = root(fobj, z, method='lm')
    steady_state = sol.x

    s = steady_state[:len(s0)]
    x = steady_state[len(s0):-n_e]
    e = steady_state[-n_e:]

    calib = OrderedDict(states=s, controls=x, shocks=e, parameters=p.copy())

    if 'auxiliary' in model.functions:
        a = model.functions['auxiliary'](s, x, p)
        calib['auxiliaries'] = a

    from dolo.compiler.misc import CalibrationDict
    return CalibrationDict(model.symbols, calib)

def residuals(model, calib=None):
    '''
    Computes the residuals associated to a calibration.

    Parameters
    ----------
    model: NumericModel

    calib: OrderedDict
        calibration dictionary (i.e. endogenous variables and parameters by
        type)

    Returns:
    --------
    OrderedDict:
        residuals vectors by equation type
    '''

    if calib is None:
        calib = model.calibration

    from collections import OrderedDict
    res = OrderedDict()

    if 'auxiliaries' not in model.symbols:

        s = calib['states']
        x = calib['controls']
        e = calib['shocks']
        p = calib['parameters']
        f = model.functions['arbitrage']
        g = model.functions['transition']

        res['transition'] = g(s, x, e, p) - s
        res['arbitrage'] = f(s, x, e, s, x, p)
    else:

        s = calib['states']
        x = calib['controls']
        y = calib['auxiliaries']
        e = calib['shocks']
        p = calib['parameters']

        f = model.functions['arbitrage']
        g = model.functions['transition']
        a = model.functions['auxiliary']

        res['transition'] = g(s, x, e, p) - s
        res['arbitrage'] = f(s, x, e, s, x, p)
        res['auxiliary'] = a(s, x, p) - y

        if 'value' in model.functions:
            rew = model.functions['value']
            v = calib['values']
            res['value'] = rew(s, x, s, x, v, p) - v

    return res
#
# def find_steady_state(model, e=None, force_states=None, constraints=None, return_jacobian=False):
#     '''n
#     Finds the steady state corresponding to exogenous shocks :math:`e`.
#
#     :param model: an "fg" model.
#     :param e: a vector with the value for the exogenous shocks.
#     :param force_values: (optional) a vector where finite values override the equilibrium conditions. For instance a vector :math:`[0,nan,nan]` would impose that the first state must be equal to 0, while the two next ones, will be determined by the model equations. This is useful, when the deterministic model has a unit root.
#     :return: a list containing a vector for the steady-states and the corresponding steady controls.
#     '''
#
#     s0 = model.calibration['states']
#     x0 = model.calibration['controls']
#     p = model.calibration['parameters']
#     z = numpy.concatenate([s0, x0])
#
#     if e is None:
#         e = numpy.zeros( len(model.symbols['shocks']) )
#     else:
#         e = e.ravel()
#
#     if constraints is not None:
#         if isinstance(constraints, (list, tuple)):
#             inds =  numpy.where( numpy.isfinite( force_values ) )[0]
#             vals = force_values[inds]
#         elif isinstance(constraints, dict):
#             inds = [model.symbols['states'].index(k) for k in force_values.keys()]
#             vals = force_values.values()
#
#     def fobj(z):
#         s = z[:len(s0)]
#         x = z[len(s0):]
#         S = model.functions['transition'](s,x,e,p)
#         r = model.functions['arbitrage'](s,x,s,x,p)
#         res = numpy.concatenate([S-s, r,  ])
#         if force_values is not None:
#             add = S[inds]-vals
#             res = numpy.concatenate([res, add])
#         return res
#
#     from trash.dolo.numeric.solver import MyJacobian
#     jac = MyJacobian(fobj)( z )
#     if return_jacobian:
#         return jac
#
#
#     rank = numpy.linalg.matrix_rank(jac)
#     if rank < len(z):
#         import warnings
#         warnings.warn("There are {} equilibrium variables to find, but the jacobian matrix is only of rank {}. The solution is indeterminate.".format(len(z),rank))
#
#     from scipy.optimize import root
#     sol = root(fobj, z, method='lm')
#     steady_state = sol.x
#
#     return [steady_state[:len(s0)], steady_state[len(s0):]]
#

if __name__ == '__main__':

    from dolo import *
    from numpy import nan

    from dolo.algos.steady_state import find_steady_state

    model = yaml_import("examples/models/open_economy.yaml")

    ss = find_steady_state(model)

    print("Steady-state variables")
    print("states: {}".format(ss[0]))
    print("controls: {}".format(ss[1]))

    jac = find_steady_state(model, return_jacobian=True)

    rank = numpy.linalg.matrix_rank(jac)

    sol2 = find_steady_state(
        model,
        force_values=[0.3, nan]
    )  # -> returns steady-state, using calibrated values as starting point
    sol3 = find_steady_state(
        model,
        force_values={'W_1': 0.3}
    )  # -> returns steady-state, using calibrated values as starting point

    print(sol2)
    print(sol3)

#    steady_state( model, e ) # -> sets exogenous values for shocks

#    steady_state( model, {'e_a':1, 'e':9}, {'k':[8,9]})

reason: import numpy,from numpy,from scipy


text:
import torch
import numpy as np
import torch.nn.functional as F
import sklearn
import sklearn.preprocessing

# https://github.com/Confusezius/Revisiting_Deep_Metric_Learning_PyTorch/blob/efddbf23ccbe267f055867b4e1c7c6693e2447c9/batchminer/rho_distance.py#L13
def distanceweightedsampling(batch, labels, lower_cutoff=0.5, upper_cutoff=1.4, contrastive_p=.2):
    """
    This methods finds all available triplets in a batch given by the classes provided in labels, and select
    triplets based on distance sampling introduced in 'Sampling Matters in Deep Embedding Learning'.

    Randomly swaps a positive sample with the negative sample for regularization.

    Args:
        batch:  np.ndarray or torch.Tensor, batch-wise embedded training samples.
        labels: np.ndarray or torch.Tensor, ground truth labels corresponding to batch.
        lower_cutoff: float, lower cutoff value for negatives that are too close to anchor embeddings. Set to literature value. They will be assigned a zero-sample probability.
        upper_cutoff: float, upper cutoff value for positives that are too far away from the anchor embeddings. Set to literature value. They will be assigned a zero-sample probability.
    Returns:
        list of sampled data tuples containing reference indices to the position IN THE BATCH.
    """
    if isinstance(labels, torch.Tensor): labels = labels.detach().cpu().numpy()
    bs = batch.shape[0]

    distances    = pdist(batch.detach()).clamp(min=lower_cutoff)

    positives, negatives = [],[]
    labels_visited = []
    anchors = []

    for i in range(bs):
        neg = labels!=labels[i]; pos = labels==labels[i]

        use_contr = np.random.choice(2, p=[1-contrastive_p, contrastive_p])
        if np.sum(pos)>1:
            anchors.append(i)
            if use_contr:
                positives.append(i)
                #Sample negatives by distance
                pos[i] = 0
                negatives.append(np.random.choice(np.where(pos)[0]))
            else:
                q_d_inv = inverse_sphere_distances(batch, distances[i], labels, labels[i])
                #Sample positives randomly
                pos[i] = 0
                positives.append(np.random.choice(np.where(pos)[0]))
                #Sample negatives by distance
                negatives.append(np.random.choice(bs,p=q_d_inv))

    sampled_triplets = [[a,p,n] for a,p,n in zip(list(range(bs)), positives, negatives)]
    return sampled_triplets

# https://github.com/Confusezius/Deep-Metric-Learning-Baselines/blob/60772745e28bc90077831bb4c9f07a233e602797/losses.py#L263
def pdist(A):
        """
        Efficient function to compute the distance matrix for a matrix A.

        Args:
            A:   Matrix/Tensor for which the distance matrix is to be computed.
            eps: float, minimal distance/clampling value to ensure no zero values.
        Returns:
            distance_matrix, clamped to ensure no zero values are passed.
        """
        prod = torch.mm(A, A.t())
        norm = prod.diag().unsqueeze(1).expand_as(prod)
        res = (norm + norm.t() - 2 * prod).clamp(min = 0)
        return res.clamp(min = 0).sqrt()

# https://github.com/Confusezius/Revisiting_Deep_Metric_Learning_PyTorch/blob/efddbf23ccbe267f055867b4e1c7c6693e2447c9/batchminer/rho_distance.py#L50
def inverse_sphere_distances(batch, dist, labels, anchor_label):
        """
        Function to utilise the distances of batch samples to compute their
        probability of occurence, and using the inverse to sample actual negatives to the resp. anchor.

        Args:
            batch:        torch.Tensor(), batch for which the sampling probabilities w.r.t to the anchor are computed. Used only to extract the shape.
            dist:         torch.Tensor(), computed distances between anchor to all batch samples.
            labels:       np.ndarray, labels for each sample for which distances were computed in dist.
            anchor_label: float, anchor label
        Returns:
            distance_matrix, clamped to ensure no zero values are passed.
        """
        bs,dim       = len(dist),batch.shape[-1]

        #negated log-distribution of distances of unit sphere in dimension <dim>
        log_q_d_inv = ((2.0 - float(dim)) * torch.log(dist) - (float(dim-3) / 2) * torch.log(1.0 - 0.25 * (dist.pow(2))))
        #Set sampling probabilities of positives to zero
        log_q_d_inv[np.where(labels==anchor_label)[0]] = 0

        q_d_inv     = torch.exp(log_q_d_inv - torch.max(log_q_d_inv)) # - max(log) for stability
        #Set sampling probabilities of positives to zero
        q_d_inv[np.where(labels==anchor_label)[0]] = 0

        ### NOTE: Cutting of values with high distances made the results slightly worse.
        # q_d_inv[np.where(dist>upper_cutoff)[0]]    = 0

        #Normalize inverted distance for probability distr.
        q_d_inv = q_d_inv/q_d_inv.sum()
        return q_d_inv.detach().cpu().numpy()

# https://github.com/Confusezius/Deep-Metric-Learning-Baselines/blob/60772745e28bc90077831bb4c9f07a233e602797/losses.py#L416
class MarginLoss(torch.nn.Module):
    def __init__(self, margin=0.2, nu=0, beta=1.2, n_classes=100, beta_constant=False, sampling_method='distance'):
        """
        Basic Margin Loss as proposed in 'Sampling Matters in Deep Embedding Learning'.

        Args:
            margin:          float, fixed triplet margin (see also TripletLoss).
            nu:              float, regularisation weight for beta. Zero by default (in literature as well).
            beta:            float, initial value for trainable class margins. Set to default literature value.
            n_classes:       int, number of target class. Required because it dictates the number of trainable class margins.
            beta_constant:   bool, set to True if betas should not be trained.
            sampling_method: str, sampling method to use to generate training triplets.
        Returns:
            Nothing!
        """
        super(MarginLoss, self).__init__()
        self.margin             = margin
        self.n_classes          = n_classes
        self.beta_constant     = beta_constant

        self.beta_val = beta
        self.beta     = beta if beta_constant else torch.nn.Parameter(torch.ones(n_classes)*beta)

        self.nu                 = nu

        self.sampling_method    = sampling_method
        self.sampler            = distanceweightedsampling

    def forward(self, batch, labels):
        """
        Args:
            batch:   torch.Tensor() [(BS x embed_dim)], batch of embeddings
            labels:  np.ndarray [(BS x 1)], for each element of the batch assigns a class [0,...,C-1]
        Returns:
            margin loss (torch.Tensor(), batch-averaged)
        """

        if isinstance(labels, torch.Tensor): labels = labels.detach().cpu().numpy()

        sampled_triplets = self.sampler(batch, labels)

        #Compute distances between anchor-positive and anchor-negative.
        d_ap, d_an = [],[]
        for triplet in sampled_triplets:
            train_triplet = {'Anchor': batch[triplet[0],:], 'Positive':batch[triplet[1],:], 'Negative':batch[triplet[2]]}

            pos_dist = ((train_triplet['Anchor']-train_triplet['Positive']).pow(2).sum()+1e-8).pow(1/2)
            neg_dist = ((train_triplet['Anchor']-train_triplet['Negative']).pow(2).sum()+1e-8).pow(1/2)

            d_ap.append(pos_dist)
            d_an.append(neg_dist)
        d_ap, d_an = torch.stack(d_ap), torch.stack(d_an)

        #Group betas together by anchor class in sampled triplets (as each beta belongs to one class).
        if self.beta_constant:
            beta = self.beta
        else:
            beta = torch.stack([self.beta[labels[triplet[0]]] for triplet in sampled_triplets]).type(torch.cuda.FloatTensor)

        #Compute actual margin postive and margin negative loss
        pos_loss = torch.nn.functional.relu(d_ap-beta+self.margin)
        neg_loss = torch.nn.functional.relu(beta-d_an+self.margin)

        #Compute normalization constant
        pair_count = torch.sum((pos_loss>0.)+(neg_loss>0.)).type(torch.cuda.FloatTensor)

        #Actual Margin Loss
        loss = torch.sum(pos_loss+neg_loss) if pair_count==0. else torch.sum(pos_loss+neg_loss)/pair_count

        #(Optional) Add regularization penalty on betas.
        if self.nu:
            beta_regularization_loss = torch.sum(beta)
            loss = loss + beta_regularisation_loss.type(torch.cuda.FloatTensor)

        return loss

# https://github.com/euwern/proxynca_pp/blob/e2fd551d90ba62f62c722782234f03fabda50320/loss.py#L52
class ProxyNCA_prob(torch.nn.Module):
    def __init__(self, nb_classes, sz_embed, scale, device):
        torch.nn.Module.__init__(self)
        self.proxies = torch.nn.Parameter(torch.randn(nb_classes, sz_embed) / 8).to(device=device)
        self.scale = scale
        self.device = device

    def forward(self, X, T):
        P = self.proxies
        #note: self.scale is equal to sqrt(1/T)
        # in the paper T = 1/9, therefore, scale = sart(1/(1/9)) = sqrt(9) = 3
        #  we need to apply sqrt because the pairwise distance is calculated as norm^2

        P = self.scale * F.normalize(P, p = 2, dim = -1)
        X = self.scale * F.normalize(X, p = 2, dim = -1)

        D = pairwise_distance(
            torch.cat(
                [X, P]
            ),
            squared = True
        )[:X.size()[0], X.size()[0]:]

        T = binarize_and_smooth_labels(
            T = T, nb_classes = len(P), smoothing_const = 0
        )

        loss = torch.sum(- T * F.log_softmax(-D, -1), -1)
        loss = loss.mean()
        return loss
# https://github.com/euwern/proxynca_pp/blob/e2fd551d90ba62f62c722782234f03fabda50320/loss.py#L9
def binarize_and_smooth_labels(T, nb_classes, smoothing_const = 0):
    device = T.device
    T = T.cpu().numpy()
    T = sklearn.preprocessing.label_binarize(
        T, classes = range(0, nb_classes)
    )
    T = T * (1 - smoothing_const)
    T[T == 0] = smoothing_const / (nb_classes - 1)
    T = torch.FloatTensor(T).to(device=device)

    return T

# https://github.com/euwern/proxynca_pp/blob/e2fd551d90ba62f62c722782234f03fabda50320/similarity.py#L6
def pairwise_distance(a, squared=False):
    """Computes the pairwise distance matrix with numerical stability."""
    pairwise_distances_squared = torch.add(
        a.pow(2).sum(dim=1, keepdim=True).expand(a.size(0), -1),
        torch.t(a).pow(2).sum(dim=0, keepdim=True).expand(a.size(0), -1)
    ) - 2 * (
        torch.mm(a, torch.t(a))
    )

    # Deal with numerical inaccuracies. Set small negatives to zero.
    pairwise_distances_squared = torch.clamp(
        pairwise_distances_squared, min=0.0
    )

    # Get the mask where the zero distances are at.
    error_mask = torch.le(pairwise_distances_squared, 0.0)
    #print(error_mask.sum())
    # Optionally take the sqrt.
    if squared:
        pairwise_distances = pairwise_distances_squared
    else:
        pairwise_distances = torch.sqrt(
            pairwise_distances_squared + error_mask.float() * 1e-16
        )

    # Undo conditionally adding 1e-16.
    pairwise_distances = torch.mul(
        pairwise_distances,
        (error_mask == False).float()
    )

    # Explicitly set diagonals to zero.
    mask_offdiagonals = 1 - torch.eye(
        *pairwise_distances.size(),
        device=pairwise_distances.device
    )
    pairwise_distances = torch.mul(pairwise_distances, mask_offdiagonals)

    return pairwise_distances

# https://github.com/Confusezius/Revisiting_Deep_Metric_Learning_PyTorch/blob/efddbf23ccbe267f055867b4e1c7c6693e2447c9/criteria/softmax.py#L12
class NormSoftmax(torch.nn.Module):
    def __init__(self, temperature, n_classes, embed_dim, loss_softmax_lr, device):
        super(NormSoftmax, self).__init__()
        self.temperature = temperature

        self.class_map = torch.nn.Parameter(torch.Tensor(n_classes, embed_dim)).to(device=device)
        stdv = 1. / np.sqrt(self.class_map.size(1))
        self.class_map.data.uniform_(-stdv, stdv)

        self.lr = loss_softmax_lr

    def forward(self, batch, labels, **kwargs):
        class_mapped_batch = torch.nn.functional.linear(batch, torch.nn.functional.normalize(self.class_map, dim=1))

        loss = torch.nn.CrossEntropyLoss()(class_mapped_batch/self.temperature, labels.to(torch.long))

        return loss


text:
import string
from itertools import product

import numpy as np
from pandas import DataFrame, MultiIndex, date_range, melt, wide_to_long
import pandas as pd

from .pandas_vb_common import setup  # noqa

class Melt(object):

    goal_time = 0.2

    def setup(self):
        self.df = DataFrame(np.random.randn(10000, 3), columns=['A', 'B', 'C'])
        self.df['id1'] = np.random.randint(0, 10, 10000)
        self.df['id2'] = np.random.randint(100, 1000, 10000)

    def time_melt_dataframe(self):
        melt(self.df, id_vars=['id1', 'id2'])

class Pivot(object):

    goal_time = 0.2

    def setup(self):
        N = 10000
        index = date_range('1/1/2000', periods=N, freq='h')
        data = {'value': np.random.randn(N * 50),
                'variable': np.arange(50).repeat(N),
                'date': np.tile(index.values, 50)}
        self.df = DataFrame(data)

    def time_reshape_pivot_time_series(self):
        self.df.pivot('date', 'variable', 'value')

class SimpleReshape(object):

    goal_time = 0.2

    def setup(self):
        arrays = [np.arange(100).repeat(100),
                  np.roll(np.tile(np.arange(100), 100), 25)]
        index = MultiIndex.from_arrays(arrays)
        self.df = DataFrame(np.random.randn(10000, 4), index=index)
        self.udf = self.df.unstack(1)

    def time_stack(self):
        self.udf.stack()

    def time_unstack(self):
        self.df.unstack(1)

class Unstack(object):

    goal_time = 0.2

    def setup(self):
        m = 100
        n = 1000

        levels = np.arange(m)
        index = MultiIndex.from_product([levels] * 2)
        columns = np.arange(n)
        values = np.arange(m * m * n).reshape(m * m, n)
        self.df = DataFrame(values, index, columns)
        self.df2 = self.df.iloc[:-1]

    def time_full_product(self):
        self.df.unstack()

    def time_without_last_row(self):
        self.df2.unstack()

class SparseIndex(object):

    goal_time = 0.2

    def setup(self):
        NUM_ROWS = 1000
        self.df = DataFrame({'A': np.random.randint(50, size=NUM_ROWS),
                             'B': np.random.randint(50, size=NUM_ROWS),
                             'C': np.random.randint(-10, 10, size=NUM_ROWS),
                             'D': np.random.randint(-10, 10, size=NUM_ROWS),
                             'E': np.random.randint(10, size=NUM_ROWS),
                             'F': np.random.randn(NUM_ROWS)})
        self.df = self.df.set_index(['A', 'B', 'C', 'D', 'E'])

    def time_unstack(self):
        self.df.unstack()

class WideToLong(object):

    goal_time = 0.2

    def setup(self):
        nyrs = 20
        nidvars = 20
        N = 5000
        self.letters = list('ABCD')
        yrvars = [l + str(num)
                  for l, num in product(self.letters, range(1, nyrs + 1))]
        columns = [str(i) for i in range(nidvars)] + yrvars
        self.df = DataFrame(np.random.randn(N, nidvars + len(yrvars)),
                            columns=columns)
        self.df['id'] = self.df.index

    def time_wide_to_long_big(self):
        wide_to_long(self.df, self.letters, i='id', j='year')

class PivotTable(object):

    goal_time = 0.2

    def setup(self):
        N = 100000
        fac1 = np.array(['A', 'B', 'C'], dtype='O')
        fac2 = np.array(['one', 'two'], dtype='O')
        ind1 = np.random.randint(0, 3, size=N)
        ind2 = np.random.randint(0, 2, size=N)
        self.df = DataFrame({'key1': fac1.take(ind1),
                             'key2': fac2.take(ind2),
                             'key3': fac2.take(ind2),
                             'value1': np.random.randn(N),
                             'value2': np.random.randn(N),
                             'value3': np.random.randn(N)})

    def time_pivot_table(self):
        self.df.pivot_table(index='key1', columns=['key2', 'key3'])

class GetDummies(object):
    goal_time = 0.2

    def setup(self):
        categories = list(string.ascii_letters[:12])
        s = pd.Series(np.random.choice(categories, size=1_000_000),
                      dtype=pd.api.types.CategoricalDtype(categories))
        self.s = s

    def time_get_dummies_1d(self):
        pd.get_dummies(self.s, sparse=False)

    def time_get_dummies_1d_sparse(self):
        pd.get_dummies(self.s, sparse=True)

max_stars_repo_name: RobMulla/pandas

max_issues_repo_name: RobMulla/pandas

max_forks_repo_name: RobMulla/pandas