text:
# coding=utf-8
"""
Module containing evaluation functions suitable for judging the performance of
a fitted LightFM model.
"""

import numpy as np

from ._lightfm_fast import (CSRMatrix,
                            calculate_auc_from_rank)

__all__ = ['precision_at_k',
           'recall_at_k',
           'auc_score',
           'reciprocal_rank']

def precision_at_k(model, test_interactions, train_interactions=None,
                   k=10, user_features=None, item_features=None,
                   preserve_rows=False, num_threads=1):
    """
    Measure the precision at k metric for a model: the fraction of known
    positives in the first k positions of the ranked list of results.
    A perfect score is 1.0.

    Parameters
    ----------

    model: LightFM instance
         the model to be evaluated
    test_interactions: np.float32 csr_matrix of shape [n_users, n_items]
         Non-zero entries representing known positives in the evaluation set.
    train_interactions: np.float32 csr_matrix of shape [n_users, n_items], optional
         Non-zero entries representing known positives in the train set. These
         will be omitted from the score calculations to avoid re-recommending
         known positives.
    k: integer, optional
         The k parameter.
    user_features: np.float32 csr_matrix of shape [n_users, n_user_features], optional
         Each row contains that user's weights over features.
    item_features: np.float32 csr_matrix of shape [n_items, n_item_features], optional
         Each row contains that item's weights over features.
    preserve_rows: boolean, optional
         When False (default), the number of rows in the output will be equal
         to the number of users with interactions in the evaluation set.
         When True, the number of rows in the output will be equal to the
         number of users.
    num_threads: int, optional
         Number of parallel computation threads to use. Should
         not be higher than the number of physical cores.

    Returns
    -------

    np.array of shape [n_users with interactions or n_users,]
         Numpy array containing precision@k scores for each user. If there are
         no interactions for a given user the returned precision will be 0.
    """

    ranks = model.predict_rank(test_interactions,
                               train_interactions=train_interactions,
                               user_features=user_features,
                               item_features=item_features,
                               num_threads=num_threads)

    ranks.data = np.less(ranks.data, k, ranks.data)

    precision = np.squeeze(np.array(ranks.sum(axis=1))) / k

    if not preserve_rows:
        precision = precision[test_interactions.getnnz(axis=1) > 0]

    return precision

def recall_at_k(model, test_interactions, train_interactions=None,
                k=10, user_features=None, item_features=None,
                preserve_rows=False, num_threads=1):
    """
    Measure the recall at k metric for a model: the number of positive items in
    the first k positions of the ranked list of results divided by the number
    of positive items in the test period. A perfect score is 1.0.

    Parameters
    ----------

    model: LightFM instance
         the model to be evaluated
    test_interactions: np.float32 csr_matrix of shape [n_users, n_items]
         Non-zero entries representing known positives in the evaluation set.
    train_interactions: np.float32 csr_matrix of shape [n_users, n_items], optional
         Non-zero entries representing known positives in the train set. These
         will be omitted from the score calculations to avoid re-recommending
         known positives.
    k: integer, optional
         The k parameter.
    user_features: np.float32 csr_matrix of shape [n_users, n_user_features], optional
         Each row contains that user's weights over features.
    item_features: np.float32 csr_matrix of shape [n_items, n_item_features], optional
         Each row contains that item's weights over features.
    preserve_rows: boolean, optional
         When False (default), the number of rows in the output will be equal
         to the number of users with interactions in the evaluation set.
         When True, the number of rows in the output will be equal to the
         number of users.
    num_threads: int, optional
         Number of parallel computation threads to use. Should
         not be higher than the number of physical cores.

    Returns
    -------

    np.array of shape [n_users with interactions or n_users,]
         Numpy array containing recall@k scores for each user. If there are no
         interactions for a given user having items in the test period, the
         returned recall will be 0.
    """

    ranks = model.predict_rank(test_interactions,
                               train_interactions=train_interactions,
                               user_features=user_features,
                               item_features=item_features,
                               num_threads=num_threads)

    ranks.data = np.less(ranks.data, k, ranks.data)

    retrieved = np.squeeze(test_interactions.getnnz(axis=1))
    hit = np.squeeze(np.array(ranks.sum(axis=1)))

    if not preserve_rows:
        hit = hit[test_interactions.getnnz(axis=1) > 0]
        retrieved = retrieved[test_interactions.getnnz(axis=1) > 0]

    return hit / retrieved

def auc_score(model, test_interactions, train_interactions=None,
              user_features=None, item_features=None,
              preserve_rows=False, num_threads=1):
    """
    Measure the ROC AUC metric for a model: the probability that a randomly
    chosen positive example has a higher score than a randomly chosen negative
    example.
    A perfect score is 1.0.

    Parameters
    ----------

    model: LightFM instance
         the model to be evaluated
    test_interactions: np.float32 csr_matrix of shape [n_users, n_items]
         Non-zero entries representing known positives in the evaluation set.
    train_interactions: np.float32 csr_matrix of shape [n_users, n_items], optional
         Non-zero entries representing known positives in the train set. These
         will be omitted from the score calculations to avoid re-recommending
         known positives.
    user_features: np.float32 csr_matrix of shape [n_users, n_user_features], optional
         Each row contains that user's weights over features.
    item_features: np.float32 csr_matrix of shape [n_items, n_item_features], optional
         Each row contains that item's weights over features.
    preserve_rows: boolean, optional
         When False (default), the number of rows in the output will be equal
         to the number of users with interactions in the evaluation set.
         When True, the number of rows in the output will be equal to the
         number of users.
    num_threads: int, optional
         Number of parallel computation threads to use. Should
         not be higher than the number of physical cores.

    Returns
    -------

    np.array of shape [n_users with interactions or n_users,]
         Numpy array containing AUC scores for each user. If there are no
         interactions for a given user the returned AUC will be 0.5.
    """

    ranks = model.predict_rank(test_interactions,
                               train_interactions=train_interactions,
                               user_features=user_features,
                               item_features=item_features,
                               num_threads=num_threads)

    assert np.all(ranks.data >= 0)

    auc = np.zeros(ranks.shape[0], dtype=np.float32)

    if train_interactions is not None:
        num_train_positives = (np.squeeze(np.array(train_interactions
                                                   .getnnz(axis=1))
                                          .astype(np.int32)))
    else:
        num_train_positives = np.zeros(test_interactions.shape[0],
                                       dtype=np.int32)

    # The second argument is modified in-place, but
    # here we don't care about the inconsistency
    # introduced into the ranks matrix.
    calculate_auc_from_rank(CSRMatrix(ranks),
                            num_train_positives,
                            ranks.data,
                            auc,
                            num_threads)

    if not preserve_rows:
        auc = auc[test_interactions.getnnz(axis=1) > 0]

    return auc

def reciprocal_rank(model, test_interactions, train_interactions=None,
                    user_features=None, item_features=None,
                    preserve_rows=False, num_threads=1):
    """
    Measure the reciprocal rank metric for a model: 1 / the rank of the highest
    ranked positive example. A perfect score is 1.0.

    Parameters
    ----------

    model: LightFM instance
         the model to be evaluated
    test_interactions: np.float32 csr_matrix of shape [n_users, n_items]
         Non-zero entries representing known positives in the evaluation set.
    train_interactions: np.float32 csr_matrix of shape [n_users, n_items], optional
         Non-zero entries representing known positives in the train set. These
         will be omitted from the score calculations to avoid re-recommending
         known positives.
    user_features: np.float32 csr_matrix of shape [n_users, n_user_features], optional
         Each row contains that user's weights over features.
    item_features: np.float32 csr_matrix of shape [n_items, n_item_features], optional
         Each row contains that item's weights over features.
    preserve_rows: boolean, optional
         When False (default), the number of rows in the output will be equal
         to the number of users with interactions in the evaluation set.
         When True, the number of rows in the output will be equal to the
         number of users.
    num_threads: int, optional
         Number of parallel computation threads to use. Should
         not be higher than the number of physical cores.

    Returns
    -------

    np.array of shape [n_users with interactions or n_users,]
         Numpy array containing reciprocal rank scores for each user.
         If there are no interactions for a given user the returned value will
         be 0.0.
    """

    ranks = model.predict_rank(test_interactions,
                               train_interactions=train_interactions,
                               user_features=user_features,
                               item_features=item_features,
                               num_threads=num_threads)

    ranks.data = 1.0 / (ranks.data + 1.0)

    ranks = np.squeeze(np.array(ranks.max(axis=1).todense()))

    if not preserve_rows:
        ranks = ranks[test_interactions.getnnz(axis=1) > 0]

    return ranks

max_stars_repo_name: TopSlugger/surender-singh

max_issues_repo_name: TopSlugger/surender-singh

max_forks_repo_name: TopSlugger/surender-singh


text:
import numpy as np
from scipy.ndimage.filters import convolve
import moviepy.editor as mpy

#### 参数和约束条件
infection_rate = 0.3
incubation_rate = 0.1

dispersion_rates  = [0, 0.07, 0.03] # for S, I, R

# 该内核会模拟人类/僵尸如何用一个位置扩散至邻近位置
dispersion_kernel = np.array([
    [0.5,  1, 0.5],
    [  1, -6,   1],
    [0.5,  1, 0.5]
    ])

france = mpy.ImageClip("france_density.png").resize(width=400)
SIR = np.zeros((3,france.h, france.w),  dtype=float)
SIR[0] = france.get_frame(0).mean(axis=2)/255

start = int(0.6*france.h), int(0.737*france.w)
SIR[1,start[0], start[1]] = 0.8 # infection in center at t=0

dt = 1.0 # 一次更新=实时1个小时
hours_per_second= 7*24 # one second in the video = one week in the model
world = {'SIR':SIR, 't':0}

##### 建模

def infection(SIR, infection_rate, incubation_rate):
    """ Computes the evolution of #Sane, #Infected, #Rampaging"""
    S,I,R = SIR
    newly_infected = infection_rate*R*S
    newly_rampaging = incubation_rate*I
    dS = - newly_infected
    dI = newly_infected - newly_rampaging
    dR = newly_rampaging
    return np.array([dS, dI, dR])

def dispersion(SIR, dispersion_kernel, dispersion_rates):
    """ Computes the dispersion (spread) of people """
    return np.array( [convolve(e, dispersion_kernel, cval=0)*r
                       for (e,r) in zip(SIR, dispersion_rates)])

def update(world):
    """ spread the epidemic for one time step """
    infect = infection(world['SIR'], infection_rate, incubation_rate)
    disperse = dispersion(world['SIR'], dispersion_kernel, dispersion_rates)
    world['SIR'] += dt*( infect + disperse)
    world['t'] += dt

# 用MoviePy制作动画
def world_to_npimage(world):
    """ Converts the world's map into a RGB image for the final video."""
    coefs = np.array([2,25,25]).reshape((3,1,1))
    accentuated_world = 255*coefs*world['SIR']
    image = accentuated_world[::-1].swapaxes(0,2).swapaxes(0,1)
    return np.minimum(255, image)

def make_frame(t):
    """ Return the frame for time t """
    while world['t'] < hours_per_second*t:
        update(world)
    return world_to_npimage(world)

animation = mpy.VideoClip(make_frame, duration=25)
animation.write_videofile('test.mp4', fps=20)

ext: py

max_stars_repo_name: XRDX/pythonExamples

max_issues_repo_name: XRDX/pythonExamples

max_forks_repo_name: XRDX/pythonExamples

reason: import numpy,from scipy


text:
import numpy as np

class Ball:
    def __init__(self, position_x, position_y, velocity_x, velocity_y, radius=5,
                 canvas_width=900, canvas_height=600):
        self.position = np.array([position_x, position_y])
        self.velocity = np.array([velocity_x, velocity_y])
        self.radius = radius
        self.canvas = np.array([canvas_width, canvas_height])

    def step(self):
        self.position += self.velocity
        self.check_collisions()

    def check_collisions(self):
        if (self.position[1] + self.radius >= self.canvas[1]) and (self.velocity[1] > 0):
            self.velocity[1] *= -1
        elif (self.position[1] - self.radius <= 0) and (self.velocity[1] < 0):
            self.velocity[1] *= -1

    def check_for_goals(self):
        if ((self.position[0] + self.radius >= self.canvas[0]) and (self.velocity[0] > 0)):
            return [1, 0]
        if ((self.position[0] - self.radius <= 0) and (self.velocity[0] < 0)):
            return [0, 1]
        return [0, 0]

    def reset_position(self):
        self.position = (self.canvas/2)

reason: import numpy


text:
"""Modified from https://github.com/CSAILVision/semantic-segmentation-pytorch"""

import os

import pandas as pd
import torch
import torch.nn as nn
import torch.nn.functional as F
from scipy.io import loadmat
from torch.nn.modules import BatchNorm2d

from inference.segmentation import resnet

# constants
ARCH_ENCODER = "resnet50dilated"
ARCH_DECODER = "ppm_deepsup"
FC_DIM = 2048
NUM_CLASS = 150
CLASS_SKY = 2
EXTENDED_CLASSES = [
    26,  # Sea
    113,  # waterfall;falls
]

# Model Builder
class ModelBuilder:
    # custom weights initialization
    @staticmethod
    def weights_init(m):
        classname = m.__class__.__name__
        if classname.find('Conv') != -1:
            nn.init.kaiming_normal_(m.weight.data)
        elif classname.find('BatchNorm') != -1:
            m.weight.data.fill_(1.)
            m.bias.data.fill_(1e-4)

    @staticmethod
    def build_encoder(arch='resnet50dilated', fc_dim=512, weights=''):
        pretrained = True if len(weights) == 0 else False
        arch = arch.lower()

        if arch == 'resnet50dilated':
            orig_resnet = resnet.__dict__['resnet50'](pretrained=pretrained)
            net_encoder = ResnetDilated(orig_resnet, dilate_scale=8)
        else:
            raise Exception('Architecture undefined!')

        # encoders are usually pretrained
        # net_encoder.apply(ModelBuilder.weights_init)
        if len(weights) > 0:
            print('Loading weights for net_encoder')
            net_encoder.load_state_dict(
                torch.load(weights, map_location=lambda storage, loc: storage), strict=False)
        return net_encoder

    @staticmethod
    def build_decoder(arch='ppm_deepsup',
                      fc_dim=512, num_class=NUM_CLASS,
                      weights='', use_softmax=False):
        arch = arch.lower()
        if arch == 'ppm_deepsup':
            net_decoder = PPMDeepsup(
                num_class=num_class,
                fc_dim=fc_dim,
                use_softmax=use_softmax)
        else:
            raise Exception('Architecture undefined!')

        net_decoder.apply(ModelBuilder.weights_init)
        if len(weights) > 0:
            print('Loading weights for net_decoder')
            net_decoder.load_state_dict(
                torch.load(weights, map_location=lambda storage, loc: storage), strict=False)
        return net_decoder

    @staticmethod
    def get_decoder(path='video_models/ade20k-resnet50dilated-ppm_deepsup/decoder_epoch_20.pth'):
        return ModelBuilder.build_decoder(arch=ARCH_DECODER, fc_dim=FC_DIM, weights=path, use_softmax=True)

    @staticmethod
    def get_encoder(path='video_models/ade20k-resnet50dilated-ppm_deepsup/encoder_epoch_20.pth'):
        return ModelBuilder.build_encoder(arch=ARCH_ENCODER, fc_dim=FC_DIM, weights=path)

def conv3x3_bn_relu(in_planes, out_planes, stride=1):
    return nn.Sequential(
        nn.Conv2d(in_planes, out_planes, kernel_size=3, stride=stride, padding=1, bias=False),
        BatchNorm2d(out_planes),
        nn.ReLU(inplace=True),
    )

class NormalizeTensor:
    def __init__(self, mean, std, inplace=False):
        """Normalize a tensor image with mean and standard deviation.
        .. note::
            This transform acts out of place by default, i.e., it does not mutates the input tensor.
        See :class:`~torchvision.transforms.Normalize` for more details.
        Args:
            tensor (Tensor): Tensor image of size (C, H, W) to be normalized.
            mean (sequence): Sequence of means for each channel.
            std (sequence): Sequence of standard deviations for each channel.
            inplace(bool,optional): Bool to make this operation inplace.
        Returns:
            Tensor: Normalized Tensor image.
        """

        self.mean = mean
        self.std = std
        self.inplace = inplace

    def __call__(self, tensor):
        if not self.inplace:
            tensor = tensor.clone()

        dtype = tensor.dtype
        mean = torch.as_tensor(self.mean, dtype=dtype, device=tensor.device)
        std = torch.as_tensor(self.std, dtype=dtype, device=tensor.device)
        tensor.sub_(mean[None, :, None, None]).div_(std[None, :, None, None])
        return tensor

class SegmentationModule(nn.Module):
    def __init__(self,
                 models_dirname,
                 use_default_normalization=True,
                 ):
        super().__init__()
        models_dirname = os.path.expandvars(models_dirname)
        encoder_path = os.path.join(models_dirname, 'encoder_epoch_20.pth')
        decoder_path = os.path.join(models_dirname, 'decoder_epoch_20.pth')
        self.encoder = ModelBuilder.get_encoder(encoder_path)
        self.decoder = ModelBuilder.get_decoder(decoder_path)
        self.use_default_normalization = use_default_normalization
        self.default_normalization = NormalizeTensor(mean=[0.485, 0.456, 0.406],
                                                     std=[0.229, 0.224, 0.225])

    @property
    def feature_maps_channels(self):
        return 256 * 2**(self.return_feature_maps_level)  # 256, 512, 1024, 2048

    def forward(self, img_data, segSize=None):
        if segSize is None:
            raise NotImplementedError("Please pass segSize param. By default: (300, 300)")

        fmaps = self.encoder(img_data, return_feature_maps=True)
        pred = self.decoder(fmaps, segSize=segSize)

        return pred

    def predict(self, tensor, imgSizes=(256, ),  # (300, 375, 450, 525, 600)
                segSize=None):
        """Entry-point for segmentation. Use this methods instead of forward
        Arguments:
            tensor {torch.Tensor} -- BCHW
        Keyword Arguments:
            imgSizes {tuple or list} -- imgSizes for segmentation input.
                default: (300, 450)
                original implementation: (300, 375, 450, 525, 600)

        """
        if segSize is None:
            segSize = tensor.shape[-2:]

        with torch.no_grad():
            if self.use_default_normalization:
                tensor = self.default_normalization(tensor)

            scores = torch.zeros(1, NUM_CLASS, segSize[0], segSize[1]).cuda()

            for img_size in imgSizes:
                img_data = F.interpolate(tensor.clone(), size=img_size)

                pred_current = self.forward(img_data, segSize=segSize)

                scores = scores + pred_current / len(imgSizes)

            return scores

    def get_edges(self, t):
        edge = torch.cuda.ByteTensor(t.size()).zero_()
        edge[:, :, :, 1:] = edge[:, :, :, 1:] | (t[:, :, :, 1:] != t[:, :, :, :-1])
        edge[:, :, :, :-1] = edge[:, :, :, :-1] | (t[:, :, :, 1:] != t[:, :, :, :-1])
        edge[:, :, 1:, :] = edge[:, :, 1:, :] | (t[:, :, 1:, :] != t[:, :, :-1, :])
        edge[:, :, :-1, :] = edge[:, :, :-1, :] | (t[:, :, 1:, :] != t[:, :, :-1, :])

        if True:
            return edge.half()
        return edge.float()

# pyramid pooling, deep supervision
class PPMDeepsup(nn.Module):
    def __init__(self, num_class=NUM_CLASS, fc_dim=4096,
                 use_softmax=False, pool_scales=(1, 2, 3, 6)):
        super().__init__()
        self.use_softmax = use_softmax

        self.ppm = []
        for scale in pool_scales:
            self.ppm.append(nn.Sequential(
                nn.AdaptiveAvgPool2d(scale),
                nn.Conv2d(fc_dim, 512, kernel_size=1, bias=False),
                BatchNorm2d(512),
                nn.ReLU(inplace=True)
            ))
        self.ppm = nn.ModuleList(self.ppm)
        self.cbr_deepsup = conv3x3_bn_relu(fc_dim // 2, fc_dim // 4, 1)

        self.conv_last = nn.Sequential(
            nn.Conv2d(fc_dim + len(pool_scales) * 512, 512,
                      kernel_size=3, padding=1, bias=False),
            BatchNorm2d(512),
            nn.ReLU(inplace=True),
            nn.Dropout2d(0.1),
            nn.Conv2d(512, num_class, kernel_size=1)
        )
        self.conv_last_deepsup = nn.Conv2d(fc_dim // 4, num_class, 1, 1, 0)
        self.dropout_deepsup = nn.Dropout2d(0.1)

    def forward(self, conv_out, segSize=None):
        conv5 = conv_out[-1]

        input_size = conv5.size()
        ppm_out = [conv5]
        for pool_scale in self.ppm:
            ppm_out.append(nn.functional.interpolate(
                pool_scale(conv5),
                (input_size[2], input_size[3]),
                mode='bilinear', align_corners=False))
        ppm_out = torch.cat(ppm_out, 1)

        x = self.conv_last(ppm_out)

        if self.use_softmax:  # is True during inference
            x = nn.functional.interpolate(
                x, size=segSize, mode='bilinear', align_corners=False)
            x = nn.functional.softmax(x, dim=1)
            return x

        # deep sup
        conv4 = conv_out[-2]
        _ = self.cbr_deepsup(conv4)
        _ = self.dropout_deepsup(_)
        _ = self.conv_last_deepsup(_)

        x = nn.functional.log_softmax(x, dim=1)
        _ = nn.functional.log_softmax(_, dim=1)

        return (x, _)

# Resnet Dilated
class ResnetDilated(nn.Module):
    def __init__(self, orig_resnet, dilate_scale=8):
        super().__init__()
        from functools import partial

        if dilate_scale == 8:
            orig_resnet.layer3.apply(
                partial(self._nostride_dilate, dilate=2))
            orig_resnet.layer4.apply(
                partial(self._nostride_dilate, dilate=4))
        elif dilate_scale == 16:
            orig_resnet.layer4.apply(
                partial(self._nostride_dilate, dilate=2))

        # take pretrained resnet, except AvgPool and FC
        self.conv1 = orig_resnet.conv1
        self.bn1 = orig_resnet.bn1
        self.relu1 = orig_resnet.relu1
        self.conv2 = orig_resnet.conv2
        self.bn2 = orig_resnet.bn2
        self.relu2 = orig_resnet.relu2
        self.conv3 = orig_resnet.conv3
        self.bn3 = orig_resnet.bn3
        self.relu3 = orig_resnet.relu3
        self.maxpool = orig_resnet.maxpool
        self.layer1 = orig_resnet.layer1
        self.layer2 = orig_resnet.layer2
        self.layer3 = orig_resnet.layer3
        self.layer4 = orig_resnet.layer4

    def _nostride_dilate(self, m, dilate):
        classname = m.__class__.__name__
        if classname.find('Conv') != -1:
            # the convolution with stride
            if m.stride == (2, 2):
                m.stride = (1, 1)
                if m.kernel_size == (3, 3):
                    m.dilation = (dilate // 2, dilate // 2)
                    m.padding = (dilate // 2, dilate // 2)
            # other convoluions
            else:
                if m.kernel_size == (3, 3):
                    m.dilation = (dilate, dilate)
                    m.padding = (dilate, dilate)

    def forward(self, x, return_feature_maps=False):
        conv_out = []

        x = self.relu1(self.bn1(self.conv1(x)))
        x = self.relu2(self.bn2(self.conv2(x)))
        x = self.relu3(self.bn3(self.conv3(x)))
        x = self.maxpool(x)

        x = self.layer1(x)
        conv_out.append(x)
        x = self.layer2(x)
        conv_out.append(x)
        x = self.layer3(x)
        conv_out.append(x)
        x = self.layer4(x)
        conv_out.append(x)

        if return_feature_maps:
            return conv_out
        return [x]

max_stars_repo_name: kuzhamuratov/deep-landscape

max_issues_repo_name: kuzhamuratov/deep-landscape

max_forks_repo_name: kuzhamuratov/deep-landscape

reason: from scipy


text:
import numpy as np

##############
### Pygame ###
##############

# Window height/width
w = 1500
h = w/2

# Where to display fitness
fitness_loc = (0.9*w, 0.05*h)

##############
### Colors ###
##############

cb = (0, 0, 0)
cdg = (252, 252, 252)
clg = (70, 70, 70)

##########################
### Physical Constants ###
##########################

# arm length
ra = int(w/10)
# ball radii
rb = int(w/100)
# current angle
o0 = -np.pi/2
# angular velocity, delta theta
do = 0
# rail length and endpoints
rdx = int(w/3)
r1 = (int(w/2 - rdx), int(h/2 - ra))
r2 = (int(w/2 + rdx), int(h/2 - ra))
# ball end coordinates (0: x, 1: y)
# use this value for centered start location
a0 = int(w/2)
# use this value for random start ocation
# a0 = np.random.choice(range(r1[0], r2[0]))
a1 = int(h/2 - ra)
b0 = a0 - int(ra * np.cos(o0))
b1 = a1 - int(ra * np.sin(o0))

##########################
### Movement Constants ###
##########################

# velocities
vd = 5
vax = 0
dv = 0
adv = 0.0035
# friction due to wall (inelastic collision)
fw = 0.45
# friction due to rail (horizonatal deceleration constant)
fr = 2
# friction at the pendulum joint (fractional angular velocity retention)
fj = 0.991
# gravitational deceleration
g = -3

###################
### Menu Params ###
###################

# See menu.py
menu_params =   {'g': (g,-10, 0),
                'fr': (fr,0,10),
                'fj': (fj, 0.8,1.1),
                'fw': (fw, 0, 1)}

defaults = (w, h, fitness_loc, cb, cdg, clg, ra, rb, o0, do, rdx, r1, r2, a0, a1, b0, b1, vd, vax, dv, adv, fw, fr, fj, g, menu_params)


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

# -------------------------------------------------------------------
# Convert float/double variables in a netCDF file to 16-bit integers
# Usage:
# float2int16.py [-h] [-3] infile outfile
#
# Convert float/double to 16-bit integers
#
# positional arguments:
#   infile      Name of input netCDF file
#   outfile     Name of output file
#
# optional arguments:
#   -h, --help  show this help message and exit
#   -3          Create netCDF-3 format instead of default netCDF-4
# -------------------------------------------------------------------

# ---------------------------------------------
# Bjørn Ådlandsvik <bjorn@imr.no>
# Institute of Marine Research, Bergen, Norway
# ---------------------------------------------

# --------
# Imports
# --------

import sys
from argparse import ArgumentParser
import collections

try:
    import numpy as np
except ImportError:
    print("ERROR: numpy is not installed")
    sys.exit(1)

try:
    from netCDF4 import Dataset
except ImportError:
    print("ERROR: netcdf4-python is not installed")
    sys.exit(1)

# -------------------------
# User settings: Configure the conversion
# -------------------------

# Variables to convert into 16-bit integers
# Format: variable = Rescale(scale_factor, add_offset)
Rescale = collections.namedtuple('Rescale', ('scale_factor', 'add_offset'))
scale_dictionary = dict(
    zeta = Rescale(0.01, 0.0),
    ubar = Rescale(0.001, 0.0),
    vbar = Rescale(0.001, 0.0),
    u    = Rescale(0.001, 0.0),
    v    = Rescale(0.001, 0.0),
    temp = Rescale(0.001, 10.0),
    salt = Rescale(0.001, 30.0),
    Akt  = Rescale(0.01, 0.0),
    Aks  = Rescale(0.01, 0.0),
)

# List of variables to not copy
dont_copy = ['omega']

# --- End user settings ---

# ----------
# Constants
# ----------

UNDEF = -32767  # 1 - 2**15

# ------------------------
# Command line arguments
# ------------------------

aparser = ArgumentParser(description="Convert float/double to 16-bit integers")

aparser.add_argument('-3', dest='format', action='store_const',
                     const='NETCDF3_CLASSIC', default='NETCDF4_CLASSIC',
                     help='Create netCDF-3 format instead of default netCDF-4')

# File names
aparser.add_argument('infile', help='Name of input netCDF file')
aparser.add_argument('outfile', help='Name of output file')

args = aparser.parse_args()

# -------------------
# Inspect input file
# -------------------

f0 = Dataset(args.infile)

# Find unlimited dimension, if any
# Classic format has at most one
unlim_dim = None
numrec = 0      # Number of records
for name in f0.dimensions:
    dim = f0.dimensions[name]
    if dim.isunlimited():
        unlim_dim = name
        numrec = len(dim)
        break

# Classify input variables
all_vars = f0.variables.keys()
for v in dont_copy:
    try:
        all_vars.remove(v)
    except ValueError:
        pass  # Don't bark if variable v is missing already

record_vars = [v for v in all_vars if unlim_dim in f0.variables[v].dimensions]
nonrec_vars = [v for v in all_vars if v not in record_vars]

# -------------------------------
# Create output file
# -------------------------------

f1 = Dataset(args.outfile, mode='w', format=args.format)

# -----------------------
# Copy global attributes
# -----------------------

for name in f0.ncattrs():
    setattr(f1, name, getattr(f0, name))

# ----------------
# Copy dimensions
# ----------------

for name, dim in f0.dimensions.items():
    if dim.isunlimited():
        f1.createDimension(name, 0)
    else:
        f1.createDimension(name, len(dim))

# ---------------------
# Variable definitions
# ---------------------
for name in f0.variables:
    v0 = f0.variables[name]
    if name in scale_dictionary:
        v1 = f1.createVariable(name, 'i2', v0.dimensions,
                               fill_value=UNDEF, zlib=True)
    else:
        if '_fillValue' in v0.ncattrs():
            f1.createVariable(name, v0.dtype, v0.dimensions,
                              fill_value=v0._fillValue)
        else:
            f1.createVariable(name, v0.dtype, v0.dimensions)

# Variable attributes
for name in f0.variables:
    v0 = f0.variables[name]
    v1 = f1.variables[name]
    for att in v0.ncattrs():
        if att != "_fillValue":
            setattr(v1, att, getattr(v0, att))
    if name in scale_dictionary:
        v1.scale_factor = scale_dictionary[name].scale_factor
        v1.add_offset = scale_dictionary[name].add_offset

# ----------------
# Non-record data
# ----------------

for name in nonrec_vars:
    v0 = f0.variables[name]
    v1 = f1.variables[name]
    v0.set_auto_maskandscale(False)
    v1.set_auto_maskandscale(False)
    if (name in scale_dictionary) and (v0.dtype != dtype('int16')):
        # Convert from float/double to int 16
        scale_factor = scale_dictionary[name].scale_factor
        offset = scale_dictionary[name].add_offset
        values = (v0[:] - offset) / scale
        values[np.abs(values) > abs(UNDEF)] = scale*UNDEF + offset
        v1[:] = np.round(values).astype('int16')
    else:  # No conversion
        v1[:] = v0[:]

# Take record variables record per record
for name in record_vars:
    v0 = f0.variables[name]
    v1 = f1.variables[name]
    v0.set_auto_maskandscale(False)
    v1.set_auto_maskandscale(False)
    if (name in scale_dictionary) and (v0.dtype != np.dtype('int16')):
        # Convert from float/double to int 16
        for rec in xrange(numrec):
            scale = scale_dictionary[name].scale_factor
            offset = scale_dictionary[name].add_offset
            values = (v0[rec, ...] - offset) / scale
            values[np.abs(values) > abs(UNDEF)] = scale*UNDEF + offset
            v1[rec, ...] = np.round(values).astype('int16')
    else:  # No conversion
        for rec in xrange(numrec):
            v1[rec, ...] = v0[rec, ...]

# ----------
# Clean up
# ----------
f1.close()

max_stars_repo_name: bjornaa/netcdf-utilites

max_issues_repo_name: bjornaa/netcdf-utilites

max_forks_repo_name: bjornaa/netcdf-utilites


text:
import numpy as np
import heapq
import random
import neuron_models as nm

neurons_per_population = 500#100
synapses_per_neuron = 400#400

delay_SP_max = 200#100
delay_SP_min = 1

weight_max = 10
weight_min = 0

fraction_negative = 0.2
ratio_inhibition_weight = 3

noise_initialization = 500
time_interval_noise = [-20,20]

time_max_sim = 3000

class LIF_Network:

    def __init__(self, neurons_per_pop = [neurons_per_population],
                 connections = [[0, 1, 400, [0, 400], [10, 10]]]):
        neuron_origin = nm.Neuron_LIF(0)
        self.count_neurons = sum(neurons_per_pop)
        self.neurons = []
        self.synapses = []
        self.spikes_heap = []
        self.neurons.append(neuron_origin)
        self.population_idx = []

        for neuron_idx in range(1,self.count_neurons +1):
            neuron = nm.Neuron_LIF(neuron_idx)
            self.neurons.append(neuron)

        curr_idx = 0
        self.population_idx.append([0,0])
        for pop_size in neurons_per_pop:
            new_idx = curr_idx + pop_size
            self.population_idx.append([curr_idx+1,new_idx])
            curr_idx = new_idx

        for con in connections:
            source = con[0]
            dest = con[1]
            synapseCount = con[2]
            time_interval = con[3]
            weight_interval = con[4]

            dest_idx_first = self.population_idx[dest][0]
            dest_idx_last = self.population_idx[dest][1]
            source_idx_first = self.population_idx[source][0]
            source_idx_last = self.population_idx[source][1]

            for dest_idx in range(dest_idx_first, dest_idx_last):
                for syn in range(synapseCount):
                    origin_idx = random.randint(source_idx_first, source_idx_last)
                    delay = random.uniform(time_interval[0], time_interval[1])
                    weight = getRandomWeight(weight_interval[0], weight_interval[1])
                    synapse = nm.Synapse( self.neurons[origin_idx], self.neurons[dest_idx], delay, weight)
                    self.synapses.append(synapse)

    def active_stdp(self, stdp_active):
        for syn in self.synapses:
            syn.stdp_active = stdp_active

    def active_learning(self, learning_binary):
        self.active_stdp(learning_binary)

    def count_active_synapses(self):
        count = 0
        for syn in self.synapses:
            if syn.active_synapse:
                count +=1
        return count

    def initial_noise(self):
        for r in range(0,noise_initialization):
            self.execute_noise(time_interval_noise)

    def execute_noise(self, time_interval):
        t_max = max(time_interval)
        t_min = min(time_interval)
        delta_t = t_min + (t_max - t_min)*np.random.rand()
        for syn in self.synapses:
            delta_t = t_min + (t_max - t_min)*np.random.rand()
            syn.stdp_noise(delta_t)

    def add_spikes(self, spike_list):
        for s in spike_list:
            heapq.heappush(self.spikes_heap,s)

    def init_activity(self, input_neuron_list):

        for neuron in input_neuron_list:
            spike = nm.generate_input_spike(neuron)
            heapq.heappush(self.spikes_heap, spike)

    def init_stimulus(self):
        spike = nm.generate_input_spike(self.neurons[0])
        heapq.heappush(self.spikes_heap, spike)

    def reset_activity(self):
        del self.spikes_heap[:]
        for n in self.neurons:
            n.reset_activity()
        for s in self.synapses:
            s.last_spike = -1

    def run(self):

        spike_count = 0
        #spikes_list_output = []
        neurons_fired = []
        firing_time = 0

        while self.spikes_heap and firing_time < time_max_sim:

            spike_count = spike_count +1
            spike = heapq.heappop(self.spikes_heap)
            #spikes_list_output.append(spike)
            execution_result = spike.execute()

            if execution_result['firing_neuron']: #meaning the receiving neuron fired

                firing_time = spike.arrival_time + 1
                firing_neuron = spike.synapse.postsynaptic_neuron
                neurons_fired.append((firing_neuron.id, firing_time))

            if execution_result['new_spikes']: #meaning new spikes were triggered
                firing_neuron = spike.synapse.postsynaptic_neuron.id
                firing_time = spike.arrival_time
                self.add_spikes(execution_result['new_spikes'])

        return neurons_fired#, 'spikes':spikes_list_output }

    def print_network(self):
        for n in self.neurons:
            n.print_neuron()

    def print_synapses(self):
        for syn in self.synapses:
            syn.print_synapse()

def get_permutation(min_int, max_int, number_elements):
    population = range(min_int, max_int)
    return random.sample(population, number_elements)

def getRandomWeight(weightMin, weightMax):
    p= random.uniform(0.,1.)
    if p<fraction_negative:
        w = random.uniform(ratio_inhibition_weight*weightMin,ratio_inhibition_weight*weightMax)
        return -w
    else:
        w = random.uniform(weightMin, weightMax)
        return w

max_stars_repo_name: pvili/SpikingTimeDependentPlasticity


text:
import os
import sys
import cv2
import numpy as np
from sklearn.cluster import MiniBatchKMeans
import math
from cPickle import dump, HIGHEST_PROTOCOL,load
from scipy.cluster.vq import vq,kmeans
import matplotlib.pyplot as plt

PRE_ALLOCATION_BUFFER=100

def getImmediateSubdirectories(dir):
	"""
		this function return the immediate subdirectory list
		eg:
			dir
				/subdirectory1
				/subdirectory2
				.
				.
		return ['subdirectory1',subdirectory2',...]
	"""

	return [name for name in os.listdir(dir) if os.path.isdir(os.path.join(dir, name))]

def getWordsInImages(imgFilesPath,useOpponent=False,disableSIFTWrite=False,Dense=False,masksPath=None,nfeatures=400):
	"""
		this function return dictonary of SIFT descriptor (simular to words in text, descriptor in image)
		they take image files path and masks path(optional)
	"""
	t=""#sys.argv[2]+"/"+"SIFT.file"
	if not(os.path.isfile(t)):
		dictWordsImages={}
		sift=cv2.SIFT(nfeatures=nfeatures)
		for imgFilepath in imgFilesPath.keys():
			if imgFilepath.endswith((".png",".jpeg",".JPG",".JPEG",".jpg",".PNG",".pgm",".PGM",".TIF",".TIFF", ".tif", ".tiff",)):
				img=cv2.imread(imgFilepath)
				width=img.shape[1]
				height=img.shape[0]
				imgGray=cv2.cvtColor(img,cv2.COLOR_BGR2GRAY)
				#imgGray=cv2.GaussianBlur(imgGray,(5,5),0)
				if masksPath!=None:
					imgMask=cv2.imread(masksPath[imgFilepath])
					imgMaskGray=cv2.cvtColor(imgMask,cv2.COLOR_BGR2GRAY)
				else:
					imgMaskGray=None

				if useOpponent and Dense:
					print "Invalid choise use either dense or useOpponent"
					sys.exit(0)
				elif useOpponent:
					d=cv2.DescriptorExtractor_create("OpponentSIFT")
					kp= sift.detect(imgGray,imgMaskGray)
					kp,des=d.compute(img,kp)
				elif Dense:
					dense=cv2.FeatureDetector_create("Dense")
					kp=dense.detect(imgGray)
					kp,des=sift.compute(imgGray,kp)
				else:
					kp, des = sift.detectAndCompute(imgGray,imgMaskGray)

				keypoints=[]
				for k in kp:
					keypoints.append(np.array([k.pt, k.size, k.angle, k.response, k.octave,k.class_id]))
				dictWordsImages[imgFilepath]=(des,keypoints,(height,width))
		if(disableSIFTWrite==False):
			file=open(t,"w")
			dump(dictWordsImages,file,protocol=HIGHEST_PROTOCOL)
	else:
		file=open(t,"r")
		dictWordsImages=load(file)
	return dictWordsImages

def dict2numpy(dict):
	"""
		this function return array from dictonary
	"""
	nkeys = len(dict)
	array = np.zeros((nkeys * PRE_ALLOCATION_BUFFER, dict[dict.keys()[0]][0].shape[1]))
	pivot = 0
	for key in dict.keys():
		(value,keypoints,shape) = dict[key]
		nelements = value.shape[0]
		while pivot + nelements > array.shape[0]:
			padding = np.zeros_like(array)
			array = np.vstack((array, padding))
		array[pivot:pivot + nelements] = value
		pivot += nelements
	array = np.resize(array, (pivot, dict[dict.keys()[0]][0].shape[1]))
	return array

def computeCentroids(arrayOfWordsImage,k,max_iter=1000,batch_size=100):
	"""
		this function compute centroid (called visual words) using k-means algorithm
	"""
	if len(sys.argv)>2:
		t=sys.argv[2]+"/"+"CENTROID.file"
	else:
		print "need to provide second command line argument"
		sys.exit(0)

	if not(os.path.isfile(t)):
		#centroids, distortion = kmeans(arrayOfWordsImage,k) used in case of scpy
		cluster=MiniBatchKMeans(init='k-means++', n_clusters=k,max_iter=max_iter,init_size=3*k,batch_size=batch_size)
		cluster.fit(arrayOfWordsImage)
		centroids = cluster.cluster_centers_
		file=open(t,"w")
		dump(centroids,file,protocol=HIGHEST_PROTOCOL)
	else:
		file=open(t,"r")
		centroids=load(file)
	return centroids

def computeHistrogram(centroids,dictWordsImages):
	"""
		this function compute the histrogram for visual words
	"""
	allwordVocabularyHistrogram={}
	for name in dictWordsImages.keys():
		(des,keypoints,shape)=dictWordsImages[name]
		if des!=None:
			code,dist=vq(des,centroids)
			#mean=(np.sum(dist)/dist.shape[0])
			#standardDevation=math.sqrt(np.sum(np.power(dist-mean,2))/dist.shape[0])
			#maxDist=standardDevation+mean
			#code=code[dist<=maxDist]
			allwordVocabularyHistrogram[name], bin_edges = np.histogram(code,bins=xrange(centroids.shape[0] + 1),normed=True)
	return allwordVocabularyHistrogram

def getCentroids():
	"""
		this function get the centroids from training data
	"""
	centroids=None
	if len(sys.argv)>4:
		t=sys.argv[4]
	else:
		print "need to provide fourth command line argument"
		sys.exit(0)
	if not(os.path.isfile(t)):
		print "Invalid path or may not be CENTROIDS.file exist USAGE:./text -d [path to testing dataset] -c [path to centroids]"
		sys.exit(0)

	else:
		try:
			file=open(t,"r")
			centroids=load(file)
			if centroids==None:
				raise Exception
		except Exception:
				import traceback
				traceback.print_exc()
	return centroids

def computeHistrogramByLevel(centroids,dictWordsImages,level=2):
	allwordVocabularyHistrogram={}
	for name in dictWordsImages:
		(des,kp,shape)=dictWordsImages[name]
		height=shape[0]
		width=shape[1]
		widthStep=width/4.0
		heightStep=height/4.0
		descriptors=des
		histogramOfLevelTwo = np.zeros((16,centroids.shape[0]))
		if descriptors!=None:
			for i in xrange(descriptors.shape[0]):
				descriptor=descriptors[i]
				keypoints=kp[i]
				x=keypoints[0][0]
				y=keypoints[0][1]
				boundaryIndex = int(x / widthStep) + int(y / heightStep) *4
				shape = descriptor.shape[0]
				descriptor = descriptor.reshape(1,shape)
				code,dist=vq(descriptor,centroids)

				histogramOfLevelTwo[boundaryIndex][code]+=1

			# level 1, based on histograms generated on level two
			histogramOfLevelOne = np.zeros((4, centroids.shape[0]))
			histogramOfLevelOne[0] = histogramOfLevelTwo[0] + histogramOfLevelTwo[1] + histogramOfLevelTwo[4] + histogramOfLevelTwo[5]
			histogramOfLevelOne[1] = histogramOfLevelTwo[2] + histogramOfLevelTwo[3] + histogramOfLevelTwo[6] + histogramOfLevelTwo[7]
			histogramOfLevelOne[2] = histogramOfLevelTwo[8] + histogramOfLevelTwo[9] + histogramOfLevelTwo[12] + histogramOfLevelTwo[13]
			histogramOfLevelOne[3] = histogramOfLevelTwo[10] + histogramOfLevelTwo[11] + histogramOfLevelTwo[14] + histogramOfLevelTwo[15]

			# level 0
			histogramOfLevelZero = histogramOfLevelOne[0] + histogramOfLevelOne[1] + histogramOfLevelOne[2] + histogramOfLevelOne[3]

			if level == 0:
				return histogramOfLevelZero
			elif level == 1:
				tempZero = histogramOfLevelZero.flatten() * 0.5
				tempOne = histogramOfLevelOne.flatten() * 0.5
				result = np.concatenate((tempZero, tempOne))
				allwordVocabularyHistrogram[name]=result
			elif level == 2:

				tempZero = histogramOfLevelZero.flatten() * 0.25
				tempOne = histogramOfLevelOne.flatten() * 0.25
				tempTwo = histogramOfLevelTwo.flatten() * 0.5
				result = np.concatenate((tempZero, tempOne, tempTwo))
				allwordVocabularyHistrogram[name]=result

			else:
				allwordVocabularyHistrogram[name]=None

	return allwordVocabularyHistrogram


text:
from __future__ import print_function
from sage.all import *

def matrix_texter(M):
    restrung=''
    b=''
    for a in str(M):
        if a=='[':
            restrung+='{'
        elif a==' ' and ( b==' ' or b=='[' or b==']'):
            restrung+=''
        elif a==' ':
            restrung+=','
        elif a==']':
            restrung+='}'
        elif a=='\n':
            restrung+=','
        else:
            restrung+=a
        b=a
    restrung='{'+restrung
    restrung=restrung+'}'
    return restrung

def vector_texter(M):
    restrung=''
    b=''
    for a in str(M):
        if a=='(':
            restrung+='{'
        elif a==')':
            restrung+='}'
        else:
            restrung+=a
        b=a
    return restrung

k_max=9
m_max=200
indix=[(k,m) for k in range(k_max+1) for m in range(m_max+1)]
for k,m in indix:
    try:
        name = 'K%d_%d' %(k,m)
        nz = load('../results/CensusKnots/'+str(k)+'T/' + name + '_invariants.sobj')
        numfield = nz[1][5]
        embedding=nz[1][7]
        shapes = nz[1][6]
        A = matrix_texter(nz[1][0])
        B = matrix_texter(nz[1][1])
        nu = nz[1][2]
        nu = vector_texter(nu)
        f = nz[1][3]
        f = vector_texter(f)
        fp = nz[1][4]
        fp = vector_texter(fp)
        textfile = open('../results/Text/CensusKnots/' + name + '.txt','w')
        print('# Manifold '+name, file=textfile)
        print(file=textfile)
        print('# Number Field', file=textfile)
        print('%s \n' %str(numfield), file=textfile)
        print('# Approximate Field Generator', file=textfile)
        print(embedding, file=textfile)
        print(file=textfile)
        print('# Shape Parameters', file=textfile)
        for shape in shapes:
            print('%s' %str(shape), file=textfile)
            print(file=textfile)
        print(file=textfile)
        print('# A Gluing Matrix', file=textfile)
        print(A, file=textfile)
        print(file=textfile)
        print('# B Gluing Matrix', file=textfile)
        print(B, file=textfile)
        print(file=textfile)
        print('# nu Gluing Vector', file=textfile)
        print(nu, file=textfile)
        print(file=textfile)
        print('# f Combinatorial flattening', file=textfile)
        print(f, file=textfile)
        print(file=textfile)
        print('# f\' Combinatorial flattening', file=textfile)
        print(fp, file=textfile)
        print(file=textfile)
        for i, invar in enumerate(nz[0]):
            print('# %d Loop Invariant' %(i+1), file=textfile)
            print(invar, file=textfile)
            print(file=textfile)
        textfile.close()
    except:
        pass

max_stars_repo_name: scottsha/compute-n-loop

max_issues_repo_name: scottsha/compute-n-loop

max_forks_repo_name: scottsha/compute-n-loop

reason: from sage


text:
"""
Copyright (c) 2019 ground0state. All rights reserved.
License: MIT License
"""
import numpy as np

class DirectionalDataAnomalyDetection():
    def __init__(self):
        self.mean_val = None

    def fit(self, X):
        """Fit the model according to the given train data.

        Parameters
        ----------
        X : array-like, shape (n_samples, n_features)
            Normal measured vectors, where n_samples is the number of samples
            and n_features is the number of features.

        Returns
        -------
        self : object
        """
        self.mean_val = X.mean(axis=0)

        return self

    def score(self, X):
        """Calculate anomaly score according to the given test data.

        Parameters
        ----------
        X: array-like, shape (n_samples, n_features)
            Error measured vectors, where n_samples is the number of samples
            and n_features is the number of features.

        Returns
        -------
        anomaly_score : array-like, shape (n_samples,)
            Anomaly score.
        """
        anomaly_score = 1 - X@self.mean_val
        return anomaly_score

if __name__ == '__main__':

    normal_data = np.loadtxt(
        "../input/normal_direction_data.csv", delimiter=",")
    error_data = np.loadtxt("../input/error_direction_data.csv", delimiter=",")

    model = DirectionalDataAnomalyDetection()
    model.fit(normal_data)
    pred = model.score(error_data)

    import matplotlib.pyplot as plt
    plt.plot(pred)
    plt.show()


text:
"""
auther: leechh
"""

import numpy as np
import tensorflow as tf
import tensorflow.keras.backend as K

def focal_loss(gamma=2., alpha=.25):
    def focal_loss_fixed(y_true, y_pred):
        pt_1 = tf.where(tf.equal(y_true, 1), y_pred, tf.ones_like(y_pred))
        pt_0 = tf.where(tf.equal(y_true, 0), y_pred, tf.zeros_like(y_pred))
        return - K.mean(alpha * K.pow(1. - pt_1, gamma) * K.log(K.epsilon()+pt_1)) \
               - K.mean((1-alpha) * K.pow(pt_0, gamma) * K.log(1. - pt_0 + K.epsilon()))
    return focal_loss_fixed

def gelu(inputs):
    return 0.5 * inputs * (1 + tf.tanh(tf.sqrt(2 / np.pi) * (inputs + 0.044715 * tf.pow(inputs, 3))))

def mish(inputs):
    return inputs * tf.math.tanh(tf.math.softplus(inputs))


text:
import os
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt

''' Read: http://pandas.pydata.org/pandas-docs/stable/api.html#api-dataframe-stats '''

def symbol_to_path(symbol, base_dir = 'data'):
	return os.path.join(base_dir, "{}.csv".format(str(symbol)))

def dates_creator(start_date, end_date):
	dates = pd.date_range(start_date, end_date)
	return dates

def get_data(symbols, dates):
	df = pd.DataFrame(index = dates)
	if 'SPY' not in symbols: # adding SPY as the main reference
		symbols.insert(0, 'SPY')
	for symbol in symbols:
		df_temp = pd.read_csv(symbol_to_path(symbol),
			index_col = 'Date',
			parse_dates = True,
			usecols = ['Date', 'Adj Close'],
			na_values = ['nan'])
		df_temp = df_temp.rename(columns = {'Adj Close': symbol})
		df = df.join(df_temp)

		if symbol == 'SPY':
			df = df.dropna(subset = ['SPY'])
	return df

def normalize_data(df):
	return df / df.iloc[0,:]

def plot(df, title):
	ax = df.plot(title = title, fontsize = 12)
	ax.set_xlabel('Date')
	ax.set_ylabel(title)
	plt.show()

def calculate_portfolio_value(normalized_df, init_investment, symbols, allocation_fraction):
        stock_investment = init_investment * allocation_fraction
        stock_values = normalized_df[symbols] * stock_investment
        portfolio_value = stock_values.sum(axis = 1)
        plot(portfolio_value, 'Portfolio value')

if __name__ == "__main__":
    start = '2013-05-01'
    end = '2013-12-31'
    symbols = ['SPY', 'AAPL', 'GOOG', 'IBM']

    dates = dates_creator(start, end)
    df = get_data(symbols, dates)
    # Daily portfolio

    normalized = normalize_data(df)

    allocation_fraction = np.array([0.5, 0.0, 0.0, 0.5])
    start_investment = 1e6 # in USD

    calculate_portfolio_value(normalized, start_investment, symbols, allocation_fraction)

max_stars_repo_name: maraboinavamshi/courses


text:
# training script

# Copyright 2020 Masao Someki
#  MIT License (https://opensource.org/licenses/MIT)

import os
import argparse

import numpy as np
import torch
from torch.utils.data import DataLoader

from net import Net
from dataset import Dataset
from net import Loss
from net import Optimizers
from writer import Logger
from utils import get_config
from decode import Decoder

np.random.seed(4)
torch.manual_seed(4)

def train(args, n_spk):
    # load config
    config = get_config(args.conf_path)

    # logger
    logger = Logger(args.log_name, 'train', 'val', 'dataset', 'decoder')

    # training device
    device = torch.device('cuda' if torch.cuda.is_available() else 'cpu')

    # hyper-set_parameters
    pad_len = 2800
    batch_len = 80
    assert pad_len % batch_len == 0

    # trainind settings and model
    net = Net(config.model, n_spk, n_cyc=2, device=device)
    net.to(device)
    iter_count = 0
    optim = Optimizers(config.optim)
    optim.set_parameters(list(net.named_parameters()))
    criteria_before = 10000
    past_model = ''

    # resume
    if args.resume is not None:
        dic = torch.load(args.resume)
        net.load_state_dict(dic['model'])
        iter_count = dic['iter_count']
        optim = dic['optim']
        criteria_before = dic['criteria']
        past_model = dic['path']

    # dataset
    datasets = {'train': Dataset(args.train_dir, args.stats_dir,
                                 logger.dataset, pad_len=pad_len,
                                 batch_len=batch_len, device=device),
                'val'  : Dataset(args.val_dir, args.stats_dir,
                                 logger.dataset, pad_len=pad_len,
                                 batch_len=batch_len, device=device)
                }

    data_loaders = {'train': DataLoader(datasets['train'],
                                        batch_size=config.train.batch_size,
                                        shuffle=True),
                    'val'  : DataLoader(datasets['val'],
                                        batch_size=config.val.batch_size,
                                        shuffle=False)}
    # logging about training data
    logger.dataset.info('number of training samples: %d' % len(datasets['train']))
    logger.dataset.info('number of validation samples: %d' % len(datasets['val']))

    # loss function
    loss_fn = Loss(device)

    # log net
    logger.train.info(net)

    # training!
    logger.train.info('Start training from iteration %d' % iter_count)

    # Hypre-parameter required to compute loss
    scale_var = torch.Tensor(datasets['train'].scaler['mcep'].scale_).to(device)

    # train
    for e in range(config.train.epoch):
        net.train()
        losses = []

        for batch in data_loaders['train']:
            # iter for batch
            iter_count += 1

            # training data to device
            inputs = {
                    'feat': torch.cat((
                                batch['uv'],
                                batch['lcf0'],
                                batch['codeap'],
                                batch['mcep']),
                            dim=-1).to(device),
                    'cv_stats': torch.cat((
                                batch['uv'],
                                batch['lcf0'],
                                batch['codeap']),
                            dim=-1).to(device),
                    'src_code': batch['src_code'].to(device),
                    'trg_code': batch['trg_code'].to(device),
                    'src_spk': batch['src_id'].to(device),
                    'trg_spk': batch['trg_id'].to(device),
                    'flen': batch['flen'].to(device)
            }

            # forward propagation
            out = net(inputs)

            # compute loss
            loss, loss_dic = loss_fn(out, inputs, scale_var)

            # backward
            loss.backward()
            optim.step()

            # log
            losses.append(loss.cpu().detach().numpy())

            if iter_count % config.train.log_every == 0:
                logger.train.info('loss at iter %d : %s' % (iter_count, loss.item()))
                logger.train.figure('train', loss_dic, iter_count)
        # log
        logger.train.info('Loss for epoch %d : %.5f' % (e, np.mean(losses)))

        # Validation
        logger.val.info('Start validation at epoch %d' % e)
        net.eval()
        losses = []
        criterias = []
        with torch.no_grad():
            for batch in data_loaders['val']:
                # to device
                inputs = {
                    'feat': torch.cat((
                            batch['uv'],
                            batch['lcf0'],
                            batch['codeap'],
                            batch['mcep']),
                        dim=-1).to(device),
                    'cv_stats': torch.cat((
                            batch['uv'],
                            batch['lcf0'],
                            batch['codeap']),
                        dim=-1).to(device),
                    'src_code': batch['src_code'].to(device),
                    'trg_code': batch['trg_code'].to(device),
                    'src_spk': batch['src_id'].to(device),
                    'trg_spk': batch['trg_id'].to(device),
                    'flen': batch['flen'].to(device)
                }

                # forward propagation
                out = net(inputs)

                # compute loss
                loss, loss_dic = loss_fn(out, inputs, scale_var)
                criteria = loss_dic['mcd/1st']

                # log
                losses.append(loss.cpu().detach().numpy())
                criterias.append(criteria)

        # log
        criteria = np.mean(criterias)
        logger.val.info('Validation loss at epoch %d: %.5f' % (e, np.mean(losses)))
        logger.val.info('Validation criteria: %.5f' % criteria)
        logger.val.figure('val', loss_dic, iter_count)

        # save model with best criteria
        if criteria < criteria_before:
            logger.val.info('Passed criteria (%f < %f), saving best model...' \
                % (criteria, criteria_before))

            # remove the existing past model.
            if not past_model == '':
                os.remove(past_model)
                logger.val.info('Found existing model at %s. Removed this file.' % past_model)

            # build dict
            save_file = os.path.join(args.model_dir, '%s.%d.pt' % (args.model_name, iter_count))
            save_dic = {
                    'model': net.state_dict(),
                    'iter_count': iter_count,
                    'optim': optim,
                    'criteria': criteria,
                    'path': save_file
                }
            torch.save(save_dic, save_file)
            criteria_before = criteria
            past_model = save_file

    # close SummaryWriter
    logger.close()

if __name__ == '__main__':
    parser = argparse.ArgumentParser()
    parser.add_argument('--train_dir', default=None, type=str,
                        help='Path to the training data.')
    parser.add_argument('--val_dir', default=None, type=str,
                        help='Path to the validation data')
    parser.add_argument('--stats_dir', default=None, type=str,
                        help='Path to the stats files.')
    parser.add_argument('--total_stats', default=None, type=str,
                        help='Path to the total_stats.h5')
    parser.add_argument('--conf_path', default=None, type=str,
                        help='Path to the config file')
    parser.add_argument('--model_dir', default=None, type=str,
                        help='Path to the directory where trained model will be saved.')
    parser.add_argument('--model_name', default=None, type=str,
                        help='Model name.')
    parser.add_argument('--decode_dir', default=None, type=str,
                        help='Path to the directory where decoded wav files will be saved')
    parser.add_argument('--log_name', default=None, type=str,
                        help='Name log file will be saved')
    parser.add_argument('--resume', default=None, type=str,
                        help='Path to the checkpoint file.')
    args = parser.parse_args()

    # number pf speakers
    n_spk = 2

    # train
    train(args, n_spk=n_spk)


text:
#
# Copyright (c) Microsoft Corporation.
# Licensed under the MIT License.
#
import math
import pytest
from typing import List

import numpy as np
import pandas as pd

from mlos.OptimizerEvaluationTools.SyntheticFunctions.Hypersphere import Hypersphere
from mlos.OptimizerEvaluationTools.ObjectiveFunctionFactory import ObjectiveFunctionFactory
from mlos.Optimizers.ParetoFrontier import ParetoFrontier
from mlos.Optimizers.OptimizationProblem import OptimizationProblem, Objective
from mlos.Spaces import ContinuousDimension, DiscreteDimension, Point, SimpleHypergrid
from mlos.Utils.KeyOrderedDict import KeyOrderedDict

class TestParetoFrontier:
    """Tests if the ParetoFrontier works."""

    def test_basic_functionality_on_2d_objective_space(self):
        """Basic sanity check. Mainly used to help us develop the API.
        """

        # Let's just create a bunch of random points, build a pareto frontier
        # and verify that the invariants hold.
        #
        parameter_space = SimpleHypergrid(
            name='params',
            dimensions=[
                ContinuousDimension(name='x1', min=0, max=10)
            ]
        )

        objective_space = SimpleHypergrid(
            name='objectives',
            dimensions=[
                ContinuousDimension(name='y1', min=0, max=10),
                ContinuousDimension(name='y2', min=0, max=10)
            ]
        )

        optimization_problem = OptimizationProblem(
            parameter_space=parameter_space,
            objective_space=objective_space,
            objectives=[
                Objective(name='y1', minimize=False),
                Objective(name='y2', minimize=False)
            ]
        )

        num_rows = 100000
        random_objectives_df = objective_space.random_dataframe(num_rows)

        # They don't match but they don't need to for this test.
        #
        random_params_df = parameter_space.random_dataframe(num_rows)

        pareto_frontier = ParetoFrontier(
            optimization_problem=optimization_problem,
            objectives_df=random_objectives_df,
            parameters_df=random_params_df
        )
        pareto_df = pareto_frontier.pareto_df

        non_pareto_index = random_objectives_df.index.difference(pareto_df.index)
        for i, row in pareto_df.iterrows():
            # Now let's make sure that no point in pareto is dominated by any non-pareto point.
            #
            assert (random_objectives_df.loc[non_pareto_index] < row).any(axis=1).sum() == len(non_pareto_index)

            # Let's also make sure that no point on the pareto is dominated by any other point there.
            #
            other_rows = pareto_df.index.difference([i])
            assert (pareto_df.loc[other_rows] > row).all(axis=1).sum() == 0

    @pytest.mark.parametrize("minimize", ["all", "none", "some"])
    @pytest.mark.parametrize("num_output_dimensions", [2, 10])
    @pytest.mark.parametrize("num_points", [100, 1000])
    def test_hyperspheres(self, minimize, num_output_dimensions, num_points):
        """Uses a hypersphere to validate that ParetoFrontier can correctly identify pareto-optimal points."""

        hypersphere_radius = 10

        objective_function_config = Point(
            implementation=Hypersphere.__name__,
            hypersphere_config=Point(
                num_objectives=num_output_dimensions,
                minimize=minimize,
                radius=hypersphere_radius
            )
        )

        objective_function = ObjectiveFunctionFactory.create_objective_function(objective_function_config=objective_function_config)
        optimization_problem = objective_function.default_optimization_problem
        random_params_df = optimization_problem.parameter_space.random_dataframe(num_points)

        # Let's randomly subsample 10% of points in random_params_df and make those points pareto optimal.
        #
        optimal_points_index = random_params_df.sample(
            frac=0.1,
            replace=False,
            axis='index'
        ).index

        random_params_df.loc[optimal_points_index, ['radius']] = hypersphere_radius
        objectives_df = objective_function.evaluate_dataframe(dataframe=random_params_df)

        # Conveniently, we can double check all of our math by invoking Pythagoras. Basically:
        #
        #   assert y0**2 + y1**2 + ... == radius**2
        #
        assert (np.power(objectives_df, 2).sum(axis=1) - np.power(random_params_df["radius"], 2) < 0.000001).all()

        # Just a few more sanity checks before we do the pareto computation.
        #
        if minimize == "all":
            assert (objectives_df <= 0).all().all()
        elif minimize == "none":
            assert (objectives_df >= 0).all().all()
        else:
            for column, minimize_column in zip(objectives_df, objective_function.minimize_mask):
                if minimize_column:
                    assert (objectives_df[column] <= 0).all()
                else:
                    assert (objectives_df[column] >= 0).all()

        pareto_frontier = ParetoFrontier(
            optimization_problem=optimization_problem,
            objectives_df=objectives_df,
            parameters_df=random_params_df
        )
        pareto_df = pareto_frontier.pareto_df

        # We know that all of the pareto efficient points must be on the frontier.
        #
        assert optimal_points_index.difference(pareto_df.index.intersection(optimal_points_index)).empty
        assert len(pareto_df.index) >= len(optimal_points_index)

        # If we flip all minimized objectives, we can assert on even more things.
        #
        for column, minimize_column in zip(objectives_df, objective_function.minimize_mask):
            if minimize_column:
                objectives_df[column] = -objectives_df[column]
                pareto_df[column] = - pareto_df[column]

        non_pareto_index = objectives_df.index.difference(pareto_df.index)
        for i, row in pareto_df.iterrows():
            # Now let's make sure that no point in pareto is dominated by any non-pareto point.
            #
            assert (objectives_df.loc[non_pareto_index] < row).any(axis=1).sum() == len(non_pareto_index)

            # Let's also make sure that no point on the pareto is dominated by any other point there.
            #
            other_rows = pareto_df.index.difference([i])
            assert (pareto_df.loc[other_rows] > row).all(axis=1).sum() == 0

    def test_repeated_values(self):
        """Validates that the algorithm does its job in the presence of repeated values.

        :return:
        """

        optimization_problem = OptimizationProblem(
            parameter_space=None,
            objective_space=SimpleHypergrid(
                name="objectives",
                dimensions=[
                    ContinuousDimension(name='y1', min=0, max=5),
                    ContinuousDimension(name='y2', min=0, max=5)
                ]
            ),
            objectives=[
                Objective(name='y1', minimize=False),
                Objective(name='y2', minimize=False)
            ]
        )

        expected_pareto_df = pd.DataFrame(
            [
                [1, 2],
                [1, 2],
                [2, 1],
                [0.5, 2],
                [1, 1],
                [2, 0.5]
            ],
            columns=['y1', 'y2']
        )

        dominated_df = pd.DataFrame(
            [
                [0.5, 0.5],
                [0.5, 1],
                [0.5, 1.5],
                [1, 0.5],
                [1.5, 0.5]
            ],
            columns=['y1', 'y2']
        )

        all_objectives_df = pd.concat([dominated_df, expected_pareto_df])
        pareto_frontier = ParetoFrontier(
            optimization_problem,
            objectives_df=all_objectives_df,
            parameters_df=pd.DataFrame(index=all_objectives_df.index)
        )
        computed_pareto_df = pareto_frontier.pareto_df
        assert computed_pareto_df.sort_values(by=['y1','y2']).equals(expected_pareto_df.sort_values(by=['y1', 'y2']))

    def test_pareto_frontier_volume_simple(self):
        """A simple sanity test on the pareto frontier volume computations.
        """

        # Let's generate a pareto frontier in 2D. ALl points lay on a line y = 1 - x
        x = np.linspace(start=0, stop=1, num=100)
        y = 1 - x
        pareto_df = pd.DataFrame({'x': x, 'y': y})
        optimization_problem = OptimizationProblem(
            parameter_space=None,
            objective_space=SimpleHypergrid(
                name='objectives',
                dimensions=[
                    ContinuousDimension(name='x', min=0, max=1),
                    ContinuousDimension(name='y', min=0, max=1)
                ]
            ),
            objectives=[Objective(name='x', minimize=False), Objective(name='y', minimize=False)]
        )
        pareto_frontier = ParetoFrontier(optimization_problem, objectives_df=pareto_df, parameters_df=pd.DataFrame(index=pareto_df.index))
        pareto_volume_estimator = pareto_frontier.approximate_pareto_volume(num_samples=1000000)
        lower_bound, upper_bound = pareto_volume_estimator.get_two_sided_confidence_interval_on_pareto_volume(alpha=0.05)
        print(lower_bound, upper_bound)
        assert 0.49 < lower_bound < upper_bound < 0.51

    @pytest.mark.parametrize("minimize", ["all", "none", "some"])
    @pytest.mark.parametrize("num_dimensions", [2, 3, 4])
    def test_pareto_frontier_volume_on_hyperspheres(self, minimize, num_dimensions):
        """Uses a known formula for the volume of the hyperspheres to validate the accuracy of the pareto frontier estimate.

        :return:
        """
        hypersphere_radius = 10
        inscribed_hypersphere_radius = 7  # For computing lower bound on volume

        # In order to validate the estimates, we must know the allowable upper and lower bounds.
        # We know that the estimate should not be higher than the volume of the n-ball (ball in n-dimensions).
        # We can also come up with a lower bound, by computing a volume of a slightly smaller ball inscribed
        # into the hypersphere. Note that the volume of an n-ball can be computed recursively, so we keep track
        # of n-ball volumes in lower dimensions.

        upper_bounds_on_sphere_volume_by_num_dimensions = {}
        lower_bounds_on_sphere_volume_by_num_dimensions = {}

        # Compute the base cases for the recursion.
        #
        upper_bounds_on_sphere_volume_by_num_dimensions[2] = np.pi * (hypersphere_radius ** 2)
        upper_bounds_on_sphere_volume_by_num_dimensions[3] = (4 / 3) * np.pi * (hypersphere_radius ** 3)

        lower_bounds_on_sphere_volume_by_num_dimensions[2] = np.pi * (inscribed_hypersphere_radius ** 2)
        lower_bounds_on_sphere_volume_by_num_dimensions[3] = (4 / 3) * np.pi * (inscribed_hypersphere_radius ** 3)

        # Compute the recursive values.
        #
        for n in range(4, num_dimensions + 1):
            upper_bounds_on_sphere_volume_by_num_dimensions[n] = upper_bounds_on_sphere_volume_by_num_dimensions[n-2] * 2 * np.pi * (hypersphere_radius ** 2) / n
            lower_bounds_on_sphere_volume_by_num_dimensions[n] = lower_bounds_on_sphere_volume_by_num_dimensions[n-2] * 2 * np.pi * (inscribed_hypersphere_radius ** 2) / n

        objective_function_config = Point(
            implementation=Hypersphere.__name__,
            hypersphere_config=Point(
                num_objectives=num_dimensions,
                minimize=minimize,
                radius=hypersphere_radius
            )
        )
        objective_function = ObjectiveFunctionFactory.create_objective_function(objective_function_config)
        parameter_space = objective_function.parameter_space

        num_points = max(4, num_dimensions)
        linspaces = []

        for dimension in parameter_space.dimensions:
            if dimension.name == 'radius':
                linspaces.append(np.array([hypersphere_radius]))
            else:
                linspaces.append(dimension.linspace(num_points))
        meshgrids = np.meshgrid(*linspaces)
        reshaped_meshgrids = [meshgrid.reshape(-1) for meshgrid in meshgrids]

        params_df = pd.DataFrame({
            dim_name: reshaped_meshgrids[i]
            for i, dim_name
            in enumerate(parameter_space.dimension_names)
        })

        objectives_df = objective_function.evaluate_dataframe(params_df)

        pareto_frontier = ParetoFrontier(
            optimization_problem=objective_function.default_optimization_problem,
            objectives_df=objectives_df,
            parameters_df=params_df
        )
        print("Num points in pareto frontier: ", len(objectives_df.index))
        assert len(pareto_frontier.pareto_df.index) == len(objectives_df.index)
        pareto_volume_estimator = pareto_frontier.approximate_pareto_volume(num_samples=1000000)
        ci_lower_bound, ci_upper_bound = pareto_volume_estimator.get_two_sided_confidence_interval_on_pareto_volume(alpha=0.05)

        lower_bound_on_pareto_volume = lower_bounds_on_sphere_volume_by_num_dimensions[num_dimensions] / (2**num_dimensions)
        upper_bound_on_pareto_volume = upper_bounds_on_sphere_volume_by_num_dimensions[num_dimensions] / (2**num_dimensions)
        print("True bounds:", lower_bound_on_pareto_volume, upper_bound_on_pareto_volume)
        print("CI bounds: ", ci_lower_bound, ci_upper_bound)
        assert lower_bound_on_pareto_volume <= ci_lower_bound <= ci_upper_bound <= upper_bound_on_pareto_volume

max_stars_repo_name: kkanellis/MLOS

max_issues_repo_name: grlap/MLOS

max_forks_repo_name: grlap/MLOS


text:
import unittest
import numpy as np
import math
import rotate3d

class TestRotate3dMisc(unittest.TestCase):
	"""A basic test suite for rotate3d helper functions"""

	def test_normalize_x(self):
		(x,y,z) = rotate3d.normalize(5,0,0)
		np.testing.assert_almost_equal(x, 1)
		np.testing.assert_almost_equal(y, 0)
		np.testing.assert_almost_equal(z, 0)

	def test_normalize_y(self):
		(x,y,z) = rotate3d.normalize(0,10,0)
		np.testing.assert_almost_equal(x, 0)
		np.testing.assert_almost_equal(y, 1)
		np.testing.assert_almost_equal(z, 0)

	def test_normalize_x(self):
		(x,y,z) = rotate3d.normalize(0,0,20)
		np.testing.assert_almost_equal(x, 0)
		np.testing.assert_almost_equal(y, 0)
		np.testing.assert_almost_equal(z, 1)

	def test_normalize_xy(self):
		(x,y,z) = rotate3d.normalize(5,5,0)
		np.testing.assert_almost_equal(x, 1/math.sqrt(2))
		np.testing.assert_almost_equal(y, 1/math.sqrt(2))
		np.testing.assert_almost_equal(z, 0)

	def test_rotation_matrix_x(self):
		matrix1 = rotate3d.axis_angle_to_rotation_matrix(1,0,0,math.radians(90))
		matrix2 = np.array([[1,0,0],[0,0,-1],[0,1,0]])
		np.testing.assert_array_almost_equal(matrix1,matrix2)

	def test_rotation_matrix_y(self):
		matrix1 = rotate3d.axis_angle_to_rotation_matrix(0,1,0,math.radians(90))
		matrix2 = np.array([[0,0,1],[0,1,0],[-1,0,0]])
		np.testing.assert_array_almost_equal(matrix1,matrix2)

class TestRotate3dObject3d(unittest.TestCase):
	"""A basic test suite for rotate3d.Object3d"""

	def setUp(self):
		v = [[1,1,-1],[1,-1,-1],[-1,-1,-1],[-1,1,-1],[0,0,1]]
		vn = []
		f = [[1,2,5],[2,3,5],[3,4,5],[4,1,5]]
		self.obj1 = rotate3d.Object3d('test_object_1', np.array(v), np.array(vn), np.array(f))
		self.obj2 = rotate3d.Object3d('test_object_2', np.array(v), np.array(vn), np.array(f))

	def test_object_init(self):
		v = [[1,1,-1],[1,-1,-1],[-1,-1,-1],[-1,1,-1],[0,0,1]]
		vn = []
		f = [[1,2,5],[2,3,5],[3,4,5],[4,1,5]]
		np.testing.assert_array_equal(self.obj1.v, v)
		np.testing.assert_array_equal(self.obj1.vn, vn)
		np.testing.assert_array_equal(self.obj1.f, f)

	def test_value_equivalence(self):
		np.testing.assert_array_equal(self.obj1.v, self.obj2.v)
		np.testing.assert_array_equal(self.obj1.vn, self.obj2.vn)
		np.testing.assert_array_equal(self.obj1.f, self.obj2.f)

	def test_rotate_and_reverse_x(self):
		self.obj1.rotate(rotate3d.axis_angle_to_rotation_matrix(1,0,0,math.radians(90)), False)
		self.obj1.rotate(rotate3d.axis_angle_to_rotation_matrix(1,0,0,math.radians(-90)), False)
		np.testing.assert_array_almost_equal(self.obj1.v, self.obj2.v)

	def test_rotate_and_reverse_y(self):
		self.obj1.rotate(rotate3d.axis_angle_to_rotation_matrix(0,1,0,math.radians(90)), False)
		self.obj1.rotate(rotate3d.axis_angle_to_rotation_matrix(0,1,0,math.radians(-90)), False)
		np.testing.assert_array_almost_equal(self.obj1.v, self.obj2.v)

	def test_rotate_and_reverse_z(self):
		self.obj1.rotate(rotate3d.axis_angle_to_rotation_matrix(0,0,1,math.radians(90)), False)
		self.obj1.rotate(rotate3d.axis_angle_to_rotation_matrix(0,0,1,math.radians(-90)), False)
		np.testing.assert_array_almost_equal(self.obj1.v, self.obj2.v)

	def test_cw_90_equals_ccw_270_x(self):
		self.obj1.rotate(rotate3d.axis_angle_to_rotation_matrix(1,0,0,math.radians(90)), False)
		self.obj2.rotate(rotate3d.axis_angle_to_rotation_matrix(1,0,0, math.radians(-270)), False)
		np.testing.assert_array_almost_equal(self.obj1.v, self.obj2.v)

	def test_cw_90_equals_ccw_270_y(self):
		self.obj1.rotate(rotate3d.axis_angle_to_rotation_matrix(0,1,0,math.radians(90)), False)
		self.obj2.rotate(rotate3d.axis_angle_to_rotation_matrix(0,1,0, math.radians(-270)), False)
		np.testing.assert_array_almost_equal(self.obj1.v, self.obj2.v)

	def test_cw_90_equals_ccw_270_z(self):
		self.obj1.rotate(rotate3d.axis_angle_to_rotation_matrix(0,0,1,math.radians(90)), False)
		self.obj2.rotate(rotate3d.axis_angle_to_rotation_matrix(0,0,1, math.radians(-270)), False)
		np.testing.assert_array_almost_equal(self.obj1.v, self.obj2.v)

	def test_colinear_axis(self):
		self.obj1.rotate(rotate3d.axis_angle_to_rotation_matrix(0,0,1,math.radians(90)), False)
		self.obj2.rotate(rotate3d.axis_angle_to_rotation_matrix(0,0,1,math.radians(90)), True)
		np.testing.assert_array_equal(self.obj1.v, self.obj2.v)

	def test_noncolinear_axis(self):
		self.obj1.rotate(rotate3d.axis_angle_to_rotation_matrix(1,0,0,math.radians(90)), True)
		self.obj2.rotate(rotate3d.axis_angle_to_rotation_matrix(1,0,0,math.radians(90)), False)
		self.assertRaises(AssertionError, np.testing.assert_array_almost_equal, self.obj1.v, self.obj2.v)

	def test_rotation_does_not_affect_faces(self):
		self.obj1.rotate(rotate3d.axis_angle_to_rotation_matrix(0,0,1,math.radians(90)), False)
		np.testing.assert_array_equal(self.obj1.f, self.obj2.f)

	def test_axis_rotations_x_y_negative_x_equals_z(self):
		self.obj1.rotate(rotate3d.axis_angle_to_rotation_matrix(1,0,0,math.radians(90)), True)
		self.obj1.rotate(rotate3d.axis_angle_to_rotation_matrix(0,1,0,math.radians(90)), True)
		self.obj1.rotate(rotate3d.axis_angle_to_rotation_matrix(1,0,0,math.radians(-90)), True)
		self.obj2.rotate(rotate3d.axis_angle_to_rotation_matrix(0,0,1,math.radians(90)), True)
		np.testing.assert_array_almost_equal(self.obj1.v, self.obj2.v)

	def test_45_degree_math_check_x(self):
		self.obj1.rotate(rotate3d.axis_angle_to_rotation_matrix(1,0,0,math.radians(45)), True)
		v = np.array([[1,0,-math.sqrt(2)],[1,-math.sqrt(2),0],[-1,-math.sqrt(2),0],[-1,0,-math.sqrt(2)],[0,1/math.sqrt(2),1/math.sqrt(2)]])
		np.testing.assert_array_almost_equal(self.obj1.v, v)

	def test_45_degree_math_check_y(self):
		self.obj1.rotate(rotate3d.axis_angle_to_rotation_matrix(0,1,0,math.radians(45)), True)
		v = np.array([[math.sqrt(2),1,0],[math.sqrt(2),-1,0],[0,-1,-math.sqrt(2)],[0,1,-math.sqrt(2)],[-1/math.sqrt(2),0,1/math.sqrt(2)]])
		np.testing.assert_array_almost_equal(self.obj1.v, v)

	def test_45_degree_math_check_z(self):
		self.obj1.rotate(rotate3d.axis_angle_to_rotation_matrix(0,0,1,math.radians(45)), True)
		v = np.array([[math.sqrt(2),0,-1],[0,-math.sqrt(2),-1],[-math.sqrt(2),0,-1],[0,math.sqrt(2),-1],[0,0,1]])
		np.testing.assert_array_almost_equal(self.obj1.v, v)

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


text:
from .abstract_learner import AbstractLearner
import numpy as np
import logging
from matplotlib import pyplot as plt
import os
import traceback
from sparse_causal_model_learner_rl.toy_datasets.dots import random_coordinates_n, image_object_positions

class DotsLearner(AbstractLearner):
    """Learn on the Dots dataset."""

    def __init__(self, *args, **kwargs):
        super().__init__(*args, **kwargs)
        self.collected = False

    def plot_images_ae(self, dataset_key='X'):
        X = self._context[dataset_key]
        idx = np.random.choice(len(X))

        img = X[idx].reshape(self.c, self.h, self.w).detach().cpu().numpy()
        pred = self.autoencoder(X[idx:idx + 1])[0].detach().cpu().numpy()

        def chw_hwc(img):
            img = np.swapaxes(img, 1, 2)
            img = np.swapaxes(img, 0, 2)
            return img

        fig = plt.figure()
        plt.subplot(1, 3, 1)
        plt.imshow(chw_hwc(img))
        plt.subplot(1, 3, 2)
        plt.imshow(chw_hwc(pred))
        plt.subplot(1, 3, 3)
        plt.imshow(np.abs(chw_hwc(img - pred)))
        return fig

    def maybe_write_artifacts(self, path_epoch, add_artifact_local):
        if (self.epochs % self.config.get('image_every', 1) == 0):
            os.makedirs(path_epoch, exist_ok=True)
            with path_epoch:
                for dataset_key in 'X_chw', 'Xtest_chw':
                    try:
                        fig = self.plot_images_ae(dataset_key=dataset_key)
                        fig.savefig(f"images_{dataset_key}.png", bbox_inches="tight")
                        artifact = path_epoch / f"images_{dataset_key}.png"
                        add_artifact_local(artifact)
                        plt.clf()
                        plt.close(fig)
                    except Exception as e:
                        logging.error(f"Images {dataset_key} error: {type(e)} {str(e)}")
                        print(traceback.format_exc())

    @property
    def _context_subclass(self):
        return {'X_chw': self.X_chw, 'n_samples': len(self.X_chw),
                'Xtest_chw': self.Xtest_chw, 'n_test': len(self.Xtest_chw)}

    def collect_steps(self):
        if self.config.get('collect_once', True):
            if self.collected:
                return
            self.collected = True

        coords_fcn = self.config.get('coords_function')

        self.X_hwc = np.array(
            [image_object_positions(positions=coords_fcn())
             for _ in range(self.config.get('n_samples_train'))])
        self.Xtest_hwc = np.array(
            [image_object_positions(positions=coords_fcn())
             for _ in range(self.config.get('n_samples_test'))])

        self.X_chw = np.swapaxes(np.swapaxes(self.X_hwc, 1, 3), 2, 3) # cwh, chw
        self.Xtest_chw = np.swapaxes(np.swapaxes(self.Xtest_hwc, 1, 3), 2, 3)

        self.h, self.w, self.c = self.X_hwc.shape[1:]
        # print("HWC", self.h, self.w, self.c)

        self.model_kwargs['input_output_shape'] = (self.h, self.w, self.c)

        print("COLLECTING DATA!!!")

    def __repr__(self):
        return f"<DotsLearner>"


text:
import pandas as pd
import numpy as np
import glob
import re
import torch
from torch import nn, optim
from torch.autograd import Variable
import torch.utils.data as data_utils
from torch.utils.data.dataset import Dataset
from gurobipy import *
import datetime
from scipy.special import softmax
from sklearn.preprocessing import OneHotEncoder
from sklearn.preprocessing import MinMaxScaler
from inspect import signature
dtype = torch.float
device = torch.device("cpu")
import logging
import warnings
## data loading
'''
separate model for each stop
'''
# npzfile = np.load("daily_stops.npz",allow_pickle=True,)
# stops = npzfile['stops_list'] # 201 length list induicating which stops are active for each day
# n_vehicles= npzfile['nr_vehicles'] #n_vehicles for each day
# weekday= npzfile['weekday'] # categorical input
# npzfile = np.load("daily_routematrix.npz", allow_pickle=True)
# opmat = npzfile['incidence_matrices'] # solutions for each day as an incidence matrix
# nxt_stops  = npzfile['next_stops']
# distance_mat = np.load("Distancematrix.npy")
# node_mat = np.load("node_category.npy")

def make_onehot(vec,num_classes):
    vec_ =  vec.reshape(len(vec),1)

    one_hot_target = (vec_ == torch.arange(num_classes).reshape(1, num_classes)).float()
    return one_hot_target

n_stops = 74
training_stops = 60
## classifier model
class costnet(nn.Module):
    def __init__(self,embedding_size,n_features=2,
                 nnodes= 74, nweekdays=7,onehot=True,**kwargs):
        super().__init__()

        if onehot:
            self.linear = nn.Linear(nnodes+nweekdays+n_features,nnodes)
        else:
            self.embeddings = nn.Embedding(nweekdays, embedding_size)
            self.linear = nn.Linear(nnodes+embedding_size+n_features,nnodes)

        self.nnodes = nnodes
        self.nweekdays = nweekdays
        self.onehot = onehot
    def forward(self,x_dist, x_features,x_week,x_mask):

        if self.onehot:
            x_emb =  make_onehot(x_week,self.nweekdays)
        else:
            x_emb = self.embeddings(x_week)

        n_batches = len(x_week)
        x = torch.cat([ x_emb, x_features, x_dist.expand(n_batches,-1)], 1)

        x = self.linear(x)
        m = nn.Softmax(dim=1)
        return m(x)* x_mask

class costnetResidual(nn.Module):
    def __init__(self, embedding_size=1,n_features=2,
                 nnodes= 74, nweekdays=7,onehot=True,**kwargs):
        super().__init__()
        if onehot:
            self.linear = nn.Linear(nweekdays+n_features,nnodes)
        else:
            self.embeddings = nn.Embedding(nweekdays, embedding_size)
            self.linear = nn.Linear(embedding_size+n_features,nnodes)
        self.nnodes = nnodes
        self.onehot = onehot
        self.nweekdays = nweekdays
    def forward(self,x_dist, x_features,x_week,x_mask):
        if self.onehot:
            x_emb =  make_onehot(x_week,self.nweekdays)
        else:
            x_emb = self.embeddings(x_week)
        x = torch.cat([ x_emb, x_features], 1)

        x = self.linear(x) - x_dist
        m = nn.Softmax(dim=1)
        return m(x)* x_mask

class VRPNet:
    def __init__(self,training_stops,
                 net=costnet,
                 embedding_size =2,nnodes=74,n_features=2, nweekdays=7,
                 optimizer=optim.Adam,epochs=20,onehot=False,**kwargs):

        self.net = net
        self.training_stops = training_stops
        self.nnodes = nnodes

        self.model = self.net(embedding_size= embedding_size,
            nnodes=nnodes,n_features=n_features, nweekdays=7,onehot=onehot)
        self.epochs = epochs

        optim_args= [k for k, v in signature(optimizer).parameters.items()]
        optim_dict = {k: kwargs.pop(k) for k in dict(kwargs) if k in optim_args}

        self.optimizer = optimizer(self.model.parameters(), **optim_dict)

    def fit(self,trgt,weekday,stops_list,n_vehicleslist,distance_mat, *arg,**kwrg):
        n_days = len(stops_list)
        x_mask = np.zeros((n_days,self.nnodes)).astype(np.float32)
        x_vehicle = np.zeros((n_days,1)).astype(np.float32)
        x_stops = np.zeros((n_days,1)).astype(np.float32)
        x_week = np.zeros(n_days).astype(np.int64)
        xdist = softmax(distance_mat[self.training_stops]).astype(np.float32)
        y_train = np.zeros((n_days,self.nnodes)).astype(np.float32)

        cnt = 0
        for d in range(n_days):
            lacti = stops_list[d]
            if self.training_stops not in lacti:
                raise Exception('Stop should be in the training days')
                # print("Stops {} Not in day {}".format(self.training_stops,
                #     d))
            else:
                x_vehicle[d] = n_vehicleslist[d]
                x_stops[d] = len(lacti)
                x_week[d] = weekday[d]
                x_mask[d][lacti] =1
                y_train[d] = trgt[d][self.training_stops]
                cnt +=1
        x_mask =  torch.from_numpy(x_mask[:cnt]).to(device)
        self.scaler = MinMaxScaler()
        x_stops = self.scaler.fit_transform(x_stops)
        x_features = np.concatenate((x_stops, x_vehicle), axis=1)
        # x_features = np.zeros((n_days,0)).astype(np.float32)

        x_features = torch.from_numpy(x_features[:cnt]).to(device)
        x_week = torch.from_numpy(x_week[:cnt]).to(device)
        xdist = torch.from_numpy(xdist).to(device)
        y_train = torch.from_numpy(y_train[:cnt]).to(device)

        criterion = nn.BCELoss() #nn.BCEWithLogitsLoss()
        for epoch in range(self.epochs):
            self.optimizer.zero_grad()
            op = self.model(xdist, x_features,x_week,x_mask)
            loss = criterion(op, y_train)
            loss.backward()
            self.optimizer.step()
            print("Epochs: {} Loss: {}".format(epoch,loss.item()))

    def predict(self,distance_mat,stops_list,weekday,n_vehicleslist,past=None, proba=True, *arg,**kwrg):
        n_days = len(stops_list)
        x_mask = np.zeros((n_days,self.nnodes)).astype(np.float32)
        x_vehicle = np.zeros((n_days,1)).astype(np.float32)
        x_stops = np.zeros((n_days,1)).astype(np.float32)
        x_week = np.zeros(n_days).astype(np.int64)
        xdist = softmax(distance_mat[self.training_stops]).astype(np.float32)
        y = np.zeros((n_days,self.nnodes)).astype(np.float32)

        cnt = 0
        active_days =[]
        for d in range(n_days):
            lacti = stops_list[d]
            if self.training_stops not in lacti:
                warnings.warn('Stop {} not in the evaluation days, \n predicted probabilities would be zero'.format(self.training_stops))

            else:
                x_vehicle[d] = n_vehicleslist[d]
                x_stops[d] = len(lacti)
                x_week[d] = weekday[d]
                x_mask[d][lacti] =1
                cnt +=1
                active_days.append(d)
        x_mask =  torch.from_numpy(x_mask[:cnt]).to(device)
        x_stops = self.scaler.transform(x_stops)
        x_features = np.concatenate((x_stops, x_vehicle), axis=1)

        x_features = torch.from_numpy(x_features[:cnt]).to(device)
        # x_vehicle = torch.from_numpy(x_vehicle[:cnt]).to(device)
        x_week = torch.from_numpy(x_week[:cnt]).to(device)
        xdist = torch.from_numpy(xdist).to(device)

        self.model.eval()
        op = self.model(xdist, x_features,x_week,x_mask)
        # m = nn.Softmax(dim=1)
        # op = m(op)

        self.model.train()
        y[active_days] = op.detach().numpy()
        if not proba:
            y = torch.argmax(y,dim=1)
        return y
    # def predict(self,stops_list,n_vehicleslist,weekday,xdist):
    #     raise Exception('Not Implemented')
class VRPresidualNet(VRPNet):
    def __init__(self,training_stops,
                 net=costnetResidual,
                 embedding_size =2,nnodes=74,n_features=2,nweekdays=7,
                 optimizer=optim.Adam,epochs=20,onehot=False, **kwargs):

        self.net = net
        self.training_stops = training_stops
        self.nnodes = nnodes

        self.model = self.net(embedding_size= embedding_size,
            nnodes=nnodes,n_features=n_features, nweekdays=7,onehot= onehot)
        self.epochs = epochs
        optim_args= [k for k, v in signature(optimizer).parameters.items()]
        optim_dict = {k: kwargs.pop(k) for k in dict(kwargs) if k in optim_args}

        self.optimizer = optimizer(self.model.parameters(), **optim_dict)