text:
from abc import abstractmethod, ABC
import numpy as np

from ptools.lipytools.little_methods import stamp
from ptools.R4C.actor import Actor
from ptools.R4C.policy_gradients.actor_critic_shared.ac_shared_graph import acs_graph
from ptools.neuralmess.nemodel import NEModel

class ACSharedModel(Actor, ABC):

    def __init__(
            self,
            observation,
            num_actions: int,
            mdict :dict,
            graph=      acs_graph,
            devices=    -1,
            verb=       1):

        self.verb = verb

        self.num_actions = num_actions
        observation_width = self.observation_vec(observation).shape[-1]

        if self.verb>0: print(f'\n*** ACSharedModel inits, observation_width: {observation_width}')

        mdict['observation_width'] = observation_width
        mdict['num_actions'] = self.num_actions
        if 'name' not in mdict: mdict['name'] = f'acs_{stamp()}'

        self.nn = NEModel(
            fwd_func=       graph,
            devices=        devices,
            save_topdir=    '_models',
            verb=           self.verb-1,
            **mdict)

        self.__upd_step = 0

    # returns vector of feats from given observation (encodes observation into input vector)
    @abstractmethod
    def observation_vec(self, observation) -> np.ndarray: pass

    # vectorization of observations batch, may be overridden with more optimal custom implementation
    def observation_vec_batch(self, observations) -> np.ndarray:
        return np.array([self.observation_vec(v) for v in observations])

    # USED for step by step actions of Actor
    def get_policy_probs(self, observation) -> np.ndarray:
        ov = self.observation_vec(observation)
        probs = self.nn.session.run(
            feed_dict=  {self.nn['observation_PH']: [ov]},
            fetches=    self.nn['action_prob'])
        return probs[0]

    def get_policy_probs_batch(self, observations) -> np.ndarray:
        ovs = self.observation_vec_batch(observations)
        probs = self.nn.session.run(
            feed_dict=  {self.nn['observation_PH']: ovs},
            fetches=    self.nn['action_prob'])
        return probs

    def get_policy_action(self, observation, sampled=False) -> int:
        probs = self.get_policy_probs(observation)
        if sampled: action = np.random.choice(self.num_actions, p=probs)
        else:       action = int(np.argmax(probs))
        return action

    def get_qvs(self, observation) -> np.ndarray:
        ov = self.observation_vec(observation)
        qvs = self.nn.session.run(
            feed_dict=  {self.nn['observation_PH']: [ov]},
            fetches=    self.nn['qvs'])
        return qvs

    def get_qvs_batch(self, observations) -> np.ndarray:
        ovs = self.observation_vec_batch(observations)
        qvss = self.nn.session.run(
            feed_dict=  {self.nn['observation_PH']: ovs},
            fetches=    self.nn['qvs'])
        return qvss

    def get_values_batch(self, observations) -> np.ndarray:
        ovs = self.observation_vec_batch(observations)
        values = self.nn.session.run(
            feed_dict=  {self.nn['observation_PH']: ovs},
            fetches=    self.nn['value'])
        return values

    def update_batch(
            self,
            observations,
            actions,
            #rewards,
            dreturns,
            next_observations,
            terminals,
            discount) -> float:

        ovs = self.observation_vec_batch(observations)

        #values_next_state = self.get_values_batch(next_observations)
        #print(values_next_state)
        #values_next_state = np.where(terminals, np.zeros_like(values_next_state), values_next_state)
        #print(values_next_state)

        #qv_labels = rewards + discount * values_next_state

        _, loss, loss_actor, loss_critic, gn, gn_avt, amax_prob, amin_prob = self.nn.session.run(
            fetches=    [
                self.nn['optimizer'],
                self.nn['loss'],
                self.nn['loss_actor'],
                self.nn['loss_critic'],
                self.nn['gg_norm'],
                self.nn['gg_avt_norm'],
                self.nn['amax_prob'],
                self.nn['amin_prob']],
            feed_dict=  {
                self.nn['observation_PH']:  ovs,
                self.nn['action_PH']:       actions,
                self.nn['qv_label_PH']:     dreturns})

        self.__upd_step += 1

        if self.verb>0:
            self.nn.log_TB(loss,        'upd/loss',         step=self.__upd_step)
            self.nn.log_TB(loss_actor,  'upd/loss_actor',   step=self.__upd_step)
            self.nn.log_TB(loss_critic, 'upd/loss_critic',  step=self.__upd_step)
            self.nn.log_TB(gn,          'upd/gn',           step=self.__upd_step)
            self.nn.log_TB(gn_avt,      'upd/gn_avt',       step=self.__upd_step)
            self.nn.log_TB(amax_prob,   'upd/amax_prob',    step=self.__upd_step)
            self.nn.log_TB(amin_prob,   'upd/amin_prob',    step=self.__upd_step)

        return loss


text:
import numpy as np
from scipy.sparse import csr_matrix

foo=np.array([[1,2,3],[4,5,6]])
print(foo)

# Sparse
foo_sparse=csr_matrix(foo)
print(foo_sparse)

# Info
print("Sparse Size: ", foo_sparse.size)
print("Sprase Shape: ", foo_sparse.shape)
print("Number of sparse dimensions: ", foo_sparse.ndim)

ext: py

max_stars_repo_name: JushBJJ/Everything

max_issues_repo_name: JushBJJ/Everything

max_forks_repo_name: JushBJJ/Everything

reason: import numpy,from scipy


text:
from enum import Enum
import numpy as np
from typing import List
import uuid

from matplotlib import pyplot as plt
from matplotlib import patches

from leap_ec.decoder import Decoder
from leap_ec.executable_rep.executable import Executable

##############################
# Class PittRulesDecoder
##############################
class PittRulesDecoder(Decoder):
    """Construct a Pitt-approach rule system (phenotype) out of a real-valued
    genome.
    """
    def __init__(self, input_space, output_space, priority_metric,
                 num_memory_registers):
        assert (input_space is not None)
        assert (output_space is not None)
        assert (num_memory_registers >= 0)
        self.input_space = input_space
        self.num_inputs = int(np.prod(input_space.shape))
        self.output_space = output_space
        self.num_outputs = int(np.prod(output_space.shape))
        self.priority_metric = priority_metric
        self.num_memory_registers = num_memory_registers

    def decode(self, genome, *args, **kwargs):
        assert (genome is not None)
        assert (len(genome) > 0)
        rule_length = self.num_inputs * 2 + \
                      self.num_outputs + self.num_memory_registers
        assert (len(genome) % rule_length == 0)
        rules = np.reshape(genome, (-1, rule_length))
        return PittRulesExecutable(self.input_space, self.output_space, rules,
                              self.priority_metric)

##############################
# Class PittRulesExecutable
##############################
class PittRulesExecutable(Executable):
    """
    An `Executable` phenotype that interprets a Pittsburgh-style ruleset and
    outputs the appropriate action.

    :param input_space: an OpenAI-gym-style space defining the inputs
    :param output_space: an OpenAI-gym-style space defining the outputs
    :param init_memory: a list of initial values for the memory registers
    :param rules: a list of rules, each of the form `[ c1 c1'  c2 c2' ... cn
        cn'  a1 ... am m1 ... mr]`
    :param priority_metric: the rule prioritization strategy used to resolve
        conflicts

    Rulesets are lists of rules.  Rules are lists of the form `[ c1 c1'  c2
    c2' ... cn cn'  a1 ... am m1 ... mr]`, where `(cx, cx')` are are the min
    and max bounds that the rule covers, `a1 .. am` are the output actions,
    and `m1 ... mr` are values to write to the memory registers.

    For example, this ruleset has two rules.  The first rule covers the
    square bounded by `(0.0, 0.6)' and `(0.0, 0.4)`, returning the output
    action `0` if the input falls within that range:

    >>> rules = [[0.0,0.6, 0.0,0.4, 0],
    ...          [0.4,1.0, 0.6,1.0, 1]]

    The input and output spaces are defined in the style of OpenAI gym.  For
    example, here's how you would set up a PittRulesExecutable with the above
    ruleset that takes two continuous input variables on `(0.0, 1.0)`,
    and outputs discrete values in `{0, 1}`:

    >>> import numpy as np
    >>> from gym import spaces
    >>> input_space = spaces.Box(low=np.array((0, 0)), high=np.array((1.0, 1.0)), dtype=np.float32)
    >>> output_space = spaces.Discrete(2)
    >>> rules = PittRulesExecutable(input_space, output_space, rules,
    ...                             priority_metric=PittRulesExecutable.PriorityMetric.RULE_ORDER)
    """

    PriorityMetric = Enum('PriorityMetric', 'RULE_ORDER GENERALITY PERIMETER')

    def __init__(self, input_space, output_space, rules, priority_metric,
                 init_mem=[]):
        assert (input_space is not None)
        assert (output_space is not None)
        assert (rules is not None)
        assert (len(rules) > 0)
        assert (priority_metric in \
                PittRulesExecutable.PriorityMetric.__members__.values())

        self.input_space = input_space
        self.num_inputs = int(np.prod(input_space.shape))
        self.output_space = output_space
        self.num_outputs = int(np.prod(output_space.shape))
        self.memory_registers = np.array(init_mem)
        self.num_memory = len(self.memory_registers)
        self.rules = rules
        self.priorities = [self.__priority(r, rule, priority_metric) for
                           (r, rule) in enumerate(rules)]

    def __priority(self, rule_order, rule, priority_metric):
        """Compute the priority value to a given rule."""
        if priority_metric == PittRulesExecutable.PriorityMetric.RULE_ORDER:
            return rule_order
        elif priority_metric == PittRulesExecutable.PriorityMetric.GENERALITY:
            pass
        elif priority_metric == PittRulesExecutable.PriorityMetric.PERIMETER:
            pass
        else:
            raise ValueError(
                'Unrecognized priority_metric "{0}".'.format(priority_metric))

    def __match_set(self, input):
        """Build the match set for a set of rules."""

        all_input = np.append(input, self.memory_registers)
        best_match_score = -1
        match_list = []  # indices of matchedrules

        for (r, rule) in enumerate(self.rules):
            match_score = 0
            for c in range(len(all_input)):  # Loop through all conditions
                # TODO Normalize this, in case the possible ranges
                # differ greatly
                low, high = (rule[c * 2], rule[c * 2 + 1])
                if low > high:
                    diff = 0  # Treat an inconsistent condition like a "wildcard": it matches antyhing
                elif all_input[c] >= low and all_input[c] <= high:
                    diff = 0  # Within range
                else:
                    diff = min(abs(low - all_input[c]), abs(high - all_input[c]))
                match_score += diff * diff  # Distance w/o sqrt

            if match_list == [] or match_score < best_match_score:
                best_match_score = match_score
                match_list = [r]
            elif match_score == best_match_score:
                match_list.append(r)
        return match_list, best_match_score

    def __fire(self, rule_index):
        winner = self.rules[rule_index]
        if self.num_memory == 0:
            return winner[-self.num_outputs:]
        else:
            self.memory_registers = winner[-self.num_memory:]
            return winner[-self.num_outputs - self.num_memory:-self.num_memory]

    def __call__(self, input_):
        """

        :param input:
        :return:

        For example, take this set of two rules:

        >>> ruleset = [[0.0,0.6, 0.0,0.5, 0],
        ...            [0.4,1.0, 0.3,1.0, 1]]

        We build an executable around it like so:

        >>> import numpy as np
        >>> from gym import spaces
        >>> input_space = spaces.Box(low=np.array((0, 0)), high=np.array((1.0, 1.0)), dtype=np.float32)
        >>> output_space = spaces.Discrete(2)
        >>> rules = PittRulesExecutable(input_space, output_space, ruleset,
        ...                             priority_metric=PittRulesExecutable.PriorityMetric.RULE_ORDER)

        It outputs `0` for inputs that are covered by only the first rule:

        >>> rules([0.1, 0.1])
        0

        >>> rules([0.5, 0.3])
        0

        It outputs `1` for inputs that are covered by only the second rule:

        >>> rules([0.9, 0.9])
        1

        >>> rules([0.5, 0.6])
        1

        If a point is covered by both rules, the first rule fires (because we set `priority_metric` to `RULE_ORDER`),
        and it outputs `0`:

        >>> rules([0.5, 0.5])
        0

        Note that if the system has more than one output, a list is returned:

        >>> ruleset = [[0.0,0.6, 0.0,0.5, 0, 1],
        ...            [0.4,1.0, 0.3,1.0, 1, 0]]
        >>> output_space = spaces.MultiBinary(2)  # A space with two binary outputs
        >>> rules = PittRulesExecutable(input_space, output_space, ruleset,
        ...                             priority_metric=PittRulesExecutable.PriorityMetric.RULE_ORDER)
        >>> rules([0.1, 0.1])
        [0, 1]

        """
        # Compute the match set
        match_list, best_match_score = self.__match_set(input_)

        # If our best-matching rules are exact matches
        if best_match_score == 0:
            # then cull the matchList based on priority.
            highest_priority_rule = min(match_list,
                                        key=lambda x: self.priorities[x])
            highest_priority_value = self.priorities[highest_priority_rule]
            match_list = [r for r in match_list if
                          self.priorities[r] <= highest_priority_value]

        # Among the rules set of rules with the highest priority, choose one
        # randomly.

        # TODO Implement other common conflict-resolutions strategies: like
        #  picking the output with the most rules

        # TODO advocating for it (i.e. vote).
        winner = np.random.choice(match_list)
        output = np.round(self.__fire(winner)).astype(int).tolist()
        assert (len(output) == self.num_outputs)
        if self.num_outputs > 1:
            return output  # Return a list of outputs if there are more than one
        else:
            return output[0]  # Return just the raw output if there is only one

##############################
# PlotPittRuleProbe class
##############################
class PlotPittRuleProbe:
    def __init__(self, context, num_inputs, num_outputs, plot_dimensions: (int, int), ax=None, xlim=(0, 1), ylim=(0, 1), modulo=1):
        assert(context is not None)
        assert(num_inputs > 0)
        assert(num_outputs > 0)
        assert(plot_dimensions is not None)
        assert(len(plot_dimensions) == 2)
        self.num_inputs = num_inputs
        self.num_outputs = num_outputs
        self.plot_dimensions = plot_dimensions
        if ax is None:
            _, ax = plt.subplots()

        ax.scatter([], [])
        self.xlim = xlim
        self.ylim = ylim
        ax.set_xlim(xlim)
        ax.set_ylim(ylim)
        self.ax = ax
        self.left, self.right = xlim
        self.bottom, self.top = ylim
        self.x = np.array([])
        self.y = np.array([])
        self.modulo = modulo
        self.context = context

    def __call__(self, population: List) -> List:
        assert(population is not None)
        assert ('leap' in self.context)
        assert ('generation' in self.context['leap'])
        step = self.context['leap']['generation']

        if step % self.modulo == 0:
            for p in reversed(self.ax.patches):
                p.remove()

            best = max(population)
            rule_length = (self.num_inputs*2 + self.num_outputs)
            if 0 != len(best.genome) % rule_length:
                raise ValueError(f"Found the wrong number of genes when trying to run {PlotPittRuleProbe.__name__}.  Are you sure these rules have {self.num_inputs} input(s) and {self.num_outputs} output(s)?")
            num_rules = int(len(best.genome)/rule_length)
            for i in range(num_rules):
                x_low = best.genome[i*rule_length + 2*self.plot_dimensions[0]]
                x_high = best.genome[i*rule_length + 2*self.plot_dimensions[0] + 1]
                y_low = best.genome[i*rule_length + 2*self.plot_dimensions[1]]
                y_high = best.genome[i*rule_length + 2*self.plot_dimensions[1] + 1]
                if x_low > x_high:
                   x_low, x_high = self.xlim
                if y_low > y_high:
                   y_low, y_high = self.ylim
                width = x_high - x_low
                height = y_high - y_low
                action = best.genome[i*rule_length + self.num_inputs]
                color = 'blue' if np.round(action) == 0 else 'red'
                self.ax.add_patch(patches.Rectangle((x_low, y_low), width, height, fill=False, color=color))

            self.ax.figure.canvas.draw()
            plt.pause(0.000001)

        return population

reason: import numpy


text:
import os

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

# from mpl_toolkits.mplot3d import Axes3D
import sklearn
from sklearn.model_selection import train_test_split
# from sklearn import svm
# from sklearn.linear_model import LinearRegression

import seaborn as sns

# TODO:Loading dataset
df = pd.read_csv('data/boston/housing.data',
             sep='\s+',
              header=None)

# Setting columns to dataset
lab_CRIM = 'Per capita crime rate by town'
lab_ZN = 'Proportion of residential land zoned for lots over 25,000 sq.ft.'
lab_INDUS = 'Proportion of non-retail business acres per town'
lab_CHAS = 'Charles River dummy variable (1 if tract bounds river; 0 otherwise)'
lab_NOX = ' nitric oxides concentration (parts per 10 million)'
lab_RM = 'Average number of rooms per dwelling'
lab_AGE = 'Proportion of owner-occupied units built prior to 1940'
lab_DIS = 'Weighted distances to five Boston employment centres'
lab_RAD = 'Index of accessibility to radial highways'
lab_TAX = 'Full-value property-tax rate per $10,000'
lab_PTRATIO = 'Pupil-teacher ratio by town'
lab_B = '1000(Bk - 0.63)^2 where Bk is the proportion of blacks by town'
lab_LSTAT = ' % lower status of the population'
lab_MEDV = 'Median value of owner-occupied homes in $1000 '

df.columns = ['CRIM', 'ZN', 'INDUS', 'CHAS', 'NOX', 'RM', 'AGE', 'DIS', 'RAD', 'TAX', 'PTRATIO', 'B', 'LSTAT', 'MEDV']

# TODO:Creating figure

prices = df.drop('MEDV', axis=1)
features = df['MEDV']

print(f'Dataset X shape: {prices.shape}')
print(f'Dataset y shape: {features.shape}')

# splitting the dataset into : train and test

x_train, x_test, y_train, y_test = train_test_split(prices, features, test_size  = 0.4 , random_state = 0)
#TODO 0.40 means 40%  of the dataset

print(x_train.shape)
print(x_test.shape)
print(y_train.shape)
print(y_test.shape)

# TODO: Print splited the dataset

print(f"x_train.shape: {x_train.shape}, y_train.shape: {y_train.shape}")
print(f"x_test.shape: {x_test.shape}, y_test.shape: {y_test.shape}")

max_stars_repo_name: cristian-nitu/final-project

max_issues_repo_name: cristian-nitu/final-project

max_forks_repo_name: cristian-nitu/final-project


text:
import os
import numpy as np
import matplotlib.pyplot as plt
from matplotlib import gridspec
from skimage.color import label2rgb
import seaborn # pylint: disable=unused-import
import cv2

import moments

def plot_image(image, fig_num, name=None, colormap='gray', write_file=None):
  """Plot image"""
  fig = plt.figure(fig_num)
  plt.imshow(image, cmap=colormap, interpolation='none')
  if name is not None:
    plt.title(name)
  plt.xticks([]), plt.yticks([])
  plt.tight_layout()
  if write_file:
    fig.savefig(write_file, bbox_inches='tight', pad_inches=0)
  fig_num += 1
  return fig_num

def plot_line_and_histogram(ordinate,
                            abscissa=None,
                            fig_num=0,
                            title='Plot and histogram',
                            label='Value',
                            num_bins=None,
                            bin_range=None,
                            normed=True,
                            marker='',
                            linestyle='.',
                            linewidth=1.0):
  """2D line plot and histogram"""

  if abscissa is None:
    abscissa = np.array(list(range(len(ordinate))))
  if num_bins is None:
    num_bins = 10
  if bin_range is None:
    bin_range = [np.min(ordinate), np.max(ordinate)]

  #fig = plt.figure(fig_num, figsize=(8, 5))
  fig = plt.figure(fig_num)
  fig.suptitle(title)
  #gs = gridspec.GridSpec(2, 1, height_ratios=[3, 1], width_ratios=[1, 1])
  gs = gridspec.GridSpec(3, 1)
  ax0 = fig.add_subplot(gs[0:2, :])
  ax0.plot(abscissa, ordinate, marker=marker, linestyle=linestyle,
           linewidth=linewidth)
  ax0.set_ylabel(label)
  ax1 = fig.add_subplot(gs[2, :])
  ax1.hist(ordinate, bins=num_bins, range=bin_range)
  ax1.set_xlabel(label)

  gs.tight_layout(fig)
  fig_num += 1

  return fig_num

def main():
  """Connected components"""
  # Configurations
  fig_num = 0
  os.makedirs('results/task02', exist_ok=True)
  # Step 1: Read image
  image_file = './images/rice.png'
  image = cv2.imread(image_file, cv2.IMREAD_GRAYSCALE)
  fig_num = plot_image(image, fig_num,
                       write_file='results/task02/original_image.png')

  # Step 2: Use morphological opening to estimate the background
  selem = cv2.getStructuringElement(cv2.MORPH_ELLIPSE, (31, 31))
  background = cv2.morphologyEx(image, cv2.MORPH_OPEN, selem)
  fig_num = plot_image(background, fig_num,
                       write_file='results/task02/background_image.png')

  # Step 3: Subtract the background image from the original image
  foreground = cv2.subtract(image, background)
  fig_num = plot_image(foreground, fig_num,
                       write_file='results/task02/foreground_image.png')

  # Step 4: Increase the image contrast
  # Here, we can implement our own, or use the built in global histogram equalization
  contrasted = cv2.equalizeHist(foreground)
  # Or we can try out the Contrast Limited Adaptive Histogram Equalization
  #clahe = cv2.createCLAHE()
  #contrasted = clahe.apply(foreground)
  fig_num = plot_image(contrasted, fig_num,
                       write_file='results/task02/histeq_image.png')
  # As is evident, none seems to work very well

  # Step 5: Threshold the image
  # The contrast adjustments seems to only make things worse (way to much
  # contrast emphasis on the background noise), so we threshold the foreground
  # image directly.
  thresh_otsu, binary_image = cv2.threshold(foreground.astype(np.uint8), 0, 255,
                                            cv2.THRESH_BINARY+cv2.THRESH_OTSU)

  print('Threshold: ', thresh_otsu)
  fig_num = plot_image(binary_image, fig_num,
                       write_file='results/task02/thresholded_image.png')

  # Step 6: Label objects in the image
  _, label_image = cv2.connectedComponents(binary_image, connectivity=4)

  # Step 7: Examine the label image
  # First, we select the grain labeled 43, and finds its contour
  single_grain_full = (label_image == 43)
  cnt_list, _ = cv2.findContours(single_grain_full.astype(np.uint8),
                                    cv2.RETR_EXTERNAL,
                                    cv2.CHAIN_APPROX_NONE)

  # As there is only one contour, we select it, and use it to compute a
  # bounding box around the grain
  grain_contour = cnt_list[0]
  y_start, x_start, dy, dx = cv2.boundingRect(grain_contour)
  grain_min_bbox = label_image[x_start:x_start+dx, y_start:y_start+dy].astype(np.uint8)

  # grain_min_bbox is the minimal bounding box around the grain, so we pad it a
  # bit before we display it.
  grain_bbox = cv2.copyMakeBorder(grain_min_bbox, 10, 10, 10, 10,
                                  cv2.BORDER_CONSTANT, value=0)
  fig_num = plot_image(grain_bbox, fig_num,
                       write_file='results/task02/grain_43.png')

  # Step 8: Display the label image
  fig_num = plot_image(label2rgb(label_image), fig_num,
                       write_file='results/task02/rgb_label_image.png')
  fig_num = plot_image(label_image, fig_num, colormap='magma',
                       write_file='results/task02/gl_label_image.png')

  # Step 9 (and 12): Measure centroid, area (and orientation) of object We compute all properties using their moments.
  # moments from interior of objects (own implementation)
  moment_dict = moments.get_moment_list(label_image)
  print(moment_dict[1])

  # Superimplose labels, centriods and orientation over this image
  fig = plt.figure(fig_num)
  plt.imshow(image, cmap='gray', interpolation='none')
  plt.xticks([]), plt.yticks([])
  plt.tight_layout()
  obj_properties = {}
  for ind, obj_moments in moment_dict.items():
    if obj_moments['m00'] > 0:
      # Area
      area = obj_moments['m00']
      theta = obj_moments['theta']
      cx = int(obj_moments['m10'] / obj_moments['m00'])
      cy = int(obj_moments['m01'] / obj_moments['m00'])
      rho = 10
      dx = rho*np.cos(theta)
      dy = rho*np.sin(theta)
      if abs(theta*180/np.pi) <= 20:
        plt.plot([cx-dx, cx, cx+dx], [cy-dy, cy, cy+dy], 'go-')
      else:
        plt.plot([cx-dx, cx, cx+dx], [cy-dy, cy, cy+dy], 'bo-')
      plt.text(cx, cy, str(ind), color='red', fontsize=14)
      # Store properties
      props = {}
      props['area'] = area
      props['cx'] = cx
      props['cy'] = cy
      props['theta'] = theta*180/np.pi
      obj_properties[ind] = props
  print(obj_properties[1])

  axes = plt.gca()
  axes.set_xlim([0, label_image.shape[0]])
  axes.set_ylim([label_image.shape[1], 0])
  fig.savefig('results/task02/marked_image.png', bbox_inches='tight', pad_inches=0)
  fig_num += 1

  areas = []
  orientations = []
  print('Object      Area Orientation Centroid')
  for ind, props in obj_properties.items():
    areas.append(props['area'])
    orientations.append(props['theta'])
    print('{0:>6} {1:>9,.0f} {2:>11,.2f} ({3:>3}, {4:>3})'.format(ind,
                                                                  props['area'],
                                                                  props['theta'],
                                                                  props['cx'],
                                                                  props['cy']))

  fig_num = plot_line_and_histogram(areas, title='Areas', marker='*',
                                    linestyle='', fig_num=fig_num, num_bins=30,
                                    normed=False, label='Area (pixels)')
  fig_num = plot_line_and_histogram(orientations, title='Orientations',
                                    marker='*', linestyle='',
                                    fig_num=fig_num, num_bins=30,
                                    normed=False, label='Angle (degrees)')

  # Index + 1 because of zero indexing in list, while the object labels are
  # enumerated from 1 (since label 0 is reserved to backround).
  print('Grain {} has the largest area of {} pixels'.format(np.argmax(areas)+1,
                                                            np.max(areas)))
  print('Mean area: {0:.2f}'.format(np.mean(areas)))

  plt.show()

if __name__ == "__main__":
  main()

max_stars_repo_name: sgzqc/wechat

max_issues_repo_name: sgzqc/wechat

max_forks_repo_name: sgzqc/wechat


text:
import torch
from scipy.stats import norm, binom_test
import numpy as np
from math import ceil
from statsmodels.stats.proportion import proportion_confint
import pandas as pd
import json
import os

# from IPython import embed

class Smooth(object):
    """A smoothed classifier g """

    # to abstain, Smooth returns this int
    ABSTAIN = -1

    def __init__(self, base_classifier: torch.nn.Module, num_classes: int, sigma: float):
        """
        :param base_classifier: maps from [batch x channel x height x width] to [batch x num_classes]
        :param num_classes:
        :param sigma: the noise level hyperparameter
        """
        self.base_classifier = base_classifier
        self.num_classes = num_classes
        self.sigma = sigma

    def certify(self, x: torch.tensor, n0: int, n: int, alpha: float, batch_size: int) -> (int, float):
        """ Monte Carlo algorithm for certifying that g's prediction around x is constant within some L2 radius.
        With probability at least 1 - alpha, the class returned by this method will equal g(x), and g's prediction will
        robust within a L2 ball of radius R around x.

        :param x: the input [channel x height x width]
        :param n0: the number of Monte Carlo samples to use for selection
        :param n: the number of Monte Carlo samples to use for estimation
        :param alpha: the failure probability
        :param batch_size: batch size to use when evaluating the base classifier
        :return: (predicted class, certified radius)
                 in the case of abstention, the class will be ABSTAIN and the radius 0.
        """
        self.base_classifier.eval()
        # draw samples of f(x+ epsilon)
        counts_selection = self._sample_noise(x, n0, batch_size)
        # use these samples to take a guess at the top class
        cAHat = counts_selection.argmax().item()
        # draw more samples of f(x + epsilon)
        counts_estimation = self._sample_noise(x, n, batch_size)
        # use these samples to estimate a lower bound on pA
        nA = counts_estimation[cAHat].item()
        pABar = self._lower_confidence_bound(nA, n, alpha)

        if pABar < 0.5:
            return Smooth.ABSTAIN, 0.0
        else:
            radius = self.sigma * norm.ppf(pABar)
            return cAHat, radius

    def predict(self, x: torch.tensor, n: int, alpha: float, batch_size: int) -> int:
        """ Monte Carlo algorithm for evaluating the prediction of g at x.  With probability at least 1 - alpha, the
        class returned by this method will equal g(x).

        This function uses the hypothesis test described in https://arxiv.org/abs/1610.03944
        for identifying the top category of a multinomial distribution.

        :param x: the input [channel x height x width]
        :param n: the number of Monte Carlo samples to use
        :param alpha: the failure probability
        :param batch_size: batch size to use when evaluating the base classifier
        :return: the predicted class, or ABSTAIN
        """
        self.base_classifier.eval()
        counts = self._sample_noise(x, n, batch_size)
        top2 = counts.argsort()[::-1][:2]
        count1 = counts[top2[0]]
        count2 = counts[top2[1]]
        if binom_test(count1, count1 + count2, p=0.5) > alpha:
            return Smooth.ABSTAIN
        else:
            return top2[0]

    def _sample_noise(self, x: torch.tensor, num: int, batch_size) -> np.ndarray:
        """ Sample the base classifier's prediction under noisy corruptions of the input x.

        :param x: the input [channel x width x height]
        :param num: number of samples to collect
        :param batch_size:
        :return: an ndarray[int] of length num_classes containing the per-class counts
        """
        with torch.no_grad():
            counts = np.zeros(self.num_classes, dtype=int)
            for _ in range(ceil(num / batch_size)):
                this_batch_size = min(batch_size, num)
                num -= this_batch_size

                batch = x.repeat((this_batch_size, 1, 1, 1))
                noise = torch.randn_like(batch, device='cuda') * self.sigma
                predictions = self.base_classifier(batch + noise).argmax(1)
                counts += self._count_arr(predictions.cpu().numpy(), self.num_classes)
            return counts

    def _count_arr(self, arr: np.ndarray, length: int) -> np.ndarray:
        counts = np.zeros(length, dtype=int)
        for idx in arr:
            counts[idx] += 1
        return counts

    def _lower_confidence_bound(self, NA: int, N: int, alpha: float) -> float:
        """ Returns a (1 - alpha) lower confidence bound on a bernoulli proportion.

        This function uses the Clopper-Pearson method.

        :param NA: the number of "successes"
        :param N: the number of total draws
        :param alpha: the confidence level
        :return: a lower bound on the binomial proportion which holds true w.p at least (1 - alpha) over the samples
        """
        return proportion_confint(NA, N, alpha=2 * alpha, method="beta")[0]

reason: import numpy,from scipy,from statsmodels


text:
import argparse
import json
import time

import jax
import numpy as np
import optax

from mesh_transformer import util
from mesh_transformer.checkpoint import read_ckpt, write_ckpt
from mesh_transformer.transformer_shard import CausalTransformer
from smart_open import open

from mesh_transformer.util import clip_by_global_norm, to_bf16

def parse_args():
    # Parse command line arguments
    parser = argparse.ArgumentParser()
    parser.add_argument("--config", type=str, default=None, help="Config file location")

    args = parser.parse_args()
    return args

if __name__ == "__main__":
    args = parse_args()
    params = json.load(open(args.config))

    cores_per_replica = params["cores_per_replica"]

    assert cores_per_replica <= 8

    bucket = params["bucket"]
    model_dir = params["model_dir"]

    params["optimizer"] = optax.chain(
        optax.scale(1),
        clip_by_global_norm(1),
        optax.scale_by_adam(),
        optax.additive_weight_decay(0),
        optax.scale(-1),
        optax.scale_by_schedule(util.gpt3_schedule(0, 1, 0, 0))
    )

    start = time.time()
    print(f"jax devices: {jax.device_count()}")
    print(f"jax runtime initialized in {time.time() - start:.06}s")

    mesh_shape = (jax.device_count() // cores_per_replica, cores_per_replica)
    devices = np.array(jax.devices()).reshape(mesh_shape)

    with open(f"gs://{bucket}/{model_dir}/meta.json", "r") as f:
        meta = json.load(f)

    ckpt_step = meta["checkpoints"][-1]
    print(f"using checkpoint {ckpt_step}")

    with jax.experimental.maps.mesh(devices, ('dp', 'mp')):
        network = CausalTransformer(params)

        start = time.time()
        network.state = read_ckpt(network.state, f"gs://{bucket}/{model_dir}/step_{ckpt_step}/", devices.shape[1])
        print(f"network loaded in {time.time() - start:.06}s")

        start = time.time()
        del network.state["opt_state"]

        network.state["params"] = to_bf16(network.state["params"])
        print(f"network converted in {time.time() - start:.06}s")

        for i in range(cores_per_replica):
            write_ckpt(network.state, f"gs://{bucket}/{model_dir}_slim/step_{ckpt_step}/", i)
            print(f"written shard {i}")

reason: import numpy,import jax


text:
import pulp
import time
from sys import maxsize
import numpy as np

def get_tracklet_cardinalities(project):
    """

    :param project: Project()
    :return: dict, {tracklet id: number of objects}
    """
    collateral_sets = get_collateral_sets(project, max_frame=5000)
    tracklets = set([t for cs in collateral_sets for t in cs])

    predecessors = generate_predecessor_map(tracklets, project)
    successors = generate_successor_map(tracklets, project)
    # print("predecessors t: {}".format(time.time() - t))

    num_animals = len(project.animals)
    t = time.time()
    areas = get_median_areas(tracklets, project)
    # TODO: medain_area should probably be median from all areas...
    median_area = np.median(list(areas.values()))
    print("area median extraction time: {}".format(time.time() - t))
    area_relaxation_coef = 0

    print("#animals {}, #tracklets: {}, #collateral sets: {}, mean area: {}".format(num_animals, len(tracklets),
                                                                                len(collateral_sets), median_area))

    cardinalities = build_ilp_and_solve(tracklets, collateral_sets, predecessors, successors,
                                                          num_animals,
                                                          median_area, areas,
                                                          print_ilp=False)  # area_relaxation_coef=area_relaxation_coef,
    return cardinalities

def fill_tracklet_cardinalites(project):
    tracklet_cardinalities = get_tracklet_cardinalities(project)
    for t_id, cardinality in list(tracklet_cardinalities.items()):
        project.chm[t_id].cardinality = cardinality
        if cardinality == 1:
            project.chm[t_id].segmentation_class = 0  # single
        elif cardinality >= 1:
            project.chm[t_id].segmentation_class = 1  # multi
        else:
            project.chm[t_id].segmentation_class = 2  # noise

def get_collateral_sets(project, max_frame=maxsize):
    max_frame = min(max_frame, project.num_frames())

    collateral_sets = []

    all_end_frames = []
    for t in project.chm.tracklet_gen():
        # if len(t) == 1:
        #     all_end_frames.append(t.end_frame())

        all_end_frames.append(t.start_frame())
        all_end_frames.append(t.end_frame())

    all_end_frames = sorted(all_end_frames)
    # all_end_frames = sorted([t.end_frame() ])

    for frame in all_end_frames:
        if frame > max_frame:
            break

        collateral_tracklets = project.chm.tracklets_in_frame(frame)
        collateral_sets.append(collateral_tracklets)

        min_end_frame = maxsize
        for tracklet in collateral_tracklets:
            min_end_frame = min(min_end_frame, tracklet.end_frame())

    return list(collateral_sets)

def get_median_areas(tracklets, project):
    import numpy as np
    areas = {}
    for t in tracklets:
        areas[t] = np.median([r.area() for r in t.r_gen(project.rm)])

    return areas

def remap_instances_to_ids(tracklets, collateral_sets, predecessors, successors, areas):
    new_tracklets = [t.id() for t in tracklets]
    new_collateral_sets = [[t.id() for t in tracklets] for tracklets in collateral_sets]
    new_predecessors = {t.id(): [t_p.id() for t_p in pred] for t, pred in list(predecessors.items())}
    new_successors = {t.id(): [t_p.id() for t_p in succ] for t, succ in list(successors.items())}
    new_areas = {t.id(): area for t, area in list(areas.items())}

    return new_tracklets, new_collateral_sets, new_predecessors, new_successors, new_areas

def relax_area_cost(cost, k, median_area, area_relaxation_coef):
    return max(0, cost - median_area * area_relaxation_coef * (k - 1))

def build_ilp(tracklets, collateral_sets, predecessors, successors, K, median_area, areas, remap_to_ids=True,
              gamma=2.0, tracklet_weights=None):
    prob = pulp.LpProblem("Cardinality classifier", pulp.LpMinimize)
    print("gamma: {}".format(gamma))

    if tracklet_weights is None:
        tracklet_weights = {t.id(): len(t) for t in tracklets}

    if remap_to_ids:
        tracklets, collateral_sets, predecessors, successors, areas = remap_instances_to_ids(tracklets,
                                                                                collateral_sets,
                                                                                predecessors,
                                                                                successors,
                                                                                areas)
    t_start = time.time()
    full_KK = list(range(0, K + 1))
    KK = list(range(1, K + 1))
    x = {}
    for t in tracklets:
        for k in full_KK:
            x[(t, k)] = pulp.LpVariable("tid:{}^{}".format(t, k), cat=pulp.LpBinary)
    print("var registration t: {}".format(time.time() - t_start))

    t_start = time.time()
    # objective function
    # prob += sum([x[(t, k)] * relax_area_cost(abs(k * median_area - areas[t]), k, median_area, area_relaxation_coef) *
    #              tracklets_len[t] for t in tracklets for k in KK]) \
    #         + sum([x[(t, 0)] * max(0, median_area - abs(median_area - areas[t])) * tracklets_len[t] * gamma for t in tracklets])  # closer to mean the worst..
    #         # + sum([x[(t, 0)] * not_used_penalty for t in tracklets])

    prob += sum([x[(t, k)] * abs(k * median_area - areas[t]) *
                 tracklet_weights[t] for t in tracklets for k in KK]) \
            + sum([x[(t, 0)] * max(0, median_area - abs(median_area - areas[t])) * tracklet_weights[t] * gamma for t in
                   tracklets])  # closer to mean the worst..
    print("objective function: {}".format(time.time() - t_start))

    t_start = time.time()
    # for each complete set
    for cs_i, cs in enumerate(collateral_sets):
        prob += sum([k * x[(t, k)] for t in cs for k in KK]) == K, "CS{}; ID perservation rule in complete sets".format(
            cs_i)

    # flow rule
    for t, pred in list(predecessors.items()):
        prob += sum([k * x[(t_, k)] for t_ in pred for k in KK]) >= \
                sum([k * x[(t, k)] for k in KK])

    for t, succ in list(successors.items()):
        prob += sum([k * x[(t_, k)] for t_ in succ for k in KK]) >= \
                sum([k * x[(t, k)] for k in KK])

    # for each tracklet - only one cardinality option is possible
    for t in tracklets:
        prob += sum([x[(t, k)] for k in full_KK]) == 1
    print("constraints: {}".format(time.time() - t_start))

    return prob, x

def generate_predecessor_map(tracklets, project):
    predecessors = {}

    for t in tracklets:
        for t_pred in t.entering_tracklets():
            if t not in predecessors:
                predecessors[t] = []

            predecessors[t].append(t_pred)

    return predecessors

def generate_successor_map(tracklets, project):
    successors = {}

    for t in tracklets:
        for t_pred in t.entering_tracklets():
            if t not in successors:
                successors[t] = []

            successors[t].append(t_pred)

    return successors

def solve(prob, x, print_ilp):
    if print_ilp:
        print(prob)

    import time

    t = time.time()
    status = prob.solve()
    print("Status:", pulp.LpStatus[prob.status])
    print("Total Cost = ", pulp.value(prob.objective))

    print("time: {}s".format(time.time() - t))
    print(status)

    cardinalities = {}
    for (id, card), var in list(x.items()):
        # id, card = str(xx)[1:-1].split(',')
        # id = int(id)
        # card = int(card[1:])
        if var.varValue > 0:
            cardinalities[id] = card
        if id not in cardinalities:
            cardinalities[id] = 0

        # print(id, card, var, var.varValue)

    # print("\n")

    # for id in sorted(cardinalities.keys()):
    #     print("{}:\t{}".format(id, cardinalities[id]))

    return cardinalities

def build_ilp_and_solve(tracklets, collateral_sets, predecesors, successors, K, median_area, areas, remap_to_ids=True,
                        print_ilp=False, gamma=2.0, tracklet_weights=None):
    """
    Build the Integral linear program and solve the cardinality problem.

    :param tracklets:
    :param collateral_sets:
    :param predecesors:
    :param successors:
    :param K:
    :param median_area:
    :param areas:
    :param remap_to_ids:
    :param print_ilp:
    :param gamma:
    :param tracklet_weights:
    :return: dict, {tracklet id: number of objects}
    """
    t = time.time()
    ilp, variable_mapping = build_ilp(tracklets, collateral_sets, predecesors, successors, K, median_area, areas,
                                      remap_to_ids,
                                      gamma, tracklet_weights=tracklet_weights)
    print("ILP construction tooks: {}s".format(time.time() - t))

    if print_ilp:
        print(ilp)

    return solve(ilp, variable_mapping, print_ilp)

def eval(tracklets, cardinalities, cardinalities_gt):
    tracklets_len = {t.id(): len(t) for t in tracklets}

    # for coef in [0, 0.1, 0.2, 0.3, 0.4, 0.5]:
    for coef in [0]:
        # cardinalities = g_cardinality_clf.build_ilp_and_solve(tracklets, collateral_sets, predecessors, num_animals,
        #                                                       median_area, areas, area_relaxation_coef=coef)
        num_mistakes = 0

        tracklets_len_sum = 0

        for id in list(cardinalities_gt.keys()):
            if cardinalities_gt[id] != cardinalities[id]:
                print("mistake, id: {}, correct cardinality: {}, estimated: {}".format(id, cardinalities_gt[id],
                                                                                cardinalities[id]))
                num_mistakes += 1
                tracklets_len_sum += tracklets_len[id]

        print("COEF: {}, #mistakes: {}, len: {}".format(coef, num_mistakes, tracklets_len_sum))
        print("------------------------------\n")

def toy_example():
    class Test():
        def __init__(self, name):
            self.id_ = name

        def id(self):
            return self.id_

    a = Test('A')
    b = Test('B')
    c = Test('C')
    d = Test('D')

    collateral_sets = [[a.id(), b.id(), c.id()], [c.id(), d.id()]]
    tracklets = set([t for cs in collateral_sets for t in cs])
    predecesors = {d.id(): [a.id(), b.id()]}
    median_area = 100
    areas = {a.id(): 110, b.id(): 85, c.id(): 400, d.id(): 180}
    K = 2

    build_ilp_and_solve(tracklets, collateral_sets, predecesors, K, median_area, areas)

if __name__ == '__main__':
    toy_example()

max_stars_repo_name: Flipajs/FERDA

max_issues_repo_name: Flipajs/FERDA

max_forks_repo_name: Flipajs/FERDA


text:
import numpy as np
import pytest
from xarray.tests import assert_allclose

from esmtools.stats import linear_slope, linregress

TIME_TYPES = ["datetime", "cftime", "float"]

@pytest.mark.parametrize("gridded", (True, False))
@pytest.mark.parametrize("time_type", TIME_TYPES)
@pytest.mark.parametrize("func", (linear_slope, linregress))
@pytest.mark.parametrize("single_arg", (True, False))
def test_linear_regression_time_da(
    gridded_da_datetime,
    gridded_da_cftime,
    gridded_da_float,
    time_type,
    gridded,
    func,
    single_arg,
):
    """Tests that linear slope can be computed on an in-memory DataArray with respect to
    time."""
    data = eval(f"gridded_da_{time_type}")()
    x = data["time"]
    y = data if gridded else data.isel(lat=0, lon=0)
    if single_arg:
        result = func(y, dim="time")
    else:
        result = func(x, y, "time")
    assert not result.isnull().any()
    if func == linear_slope:
        expected_shape = (3, 3) if gridded else ()
    elif func == linregress:
        expected_shape = (3, 3, 5) if gridded else (5,)
    assert result.shape == expected_shape

@pytest.mark.parametrize("gridded", (True, False))
@pytest.mark.parametrize("time_type", TIME_TYPES)
@pytest.mark.parametrize("func", (linear_slope, linregress))
@pytest.mark.parametrize("single_arg", (True, False))
def test_linear_regression_singular_data_da(
    gridded_da_datetime,
    gridded_da_cftime,
    gridded_da_float,
    time_type,
    gridded,
    func,
    single_arg,
):
    """Tests that ``linear_slope`` works between a grid or time series of data and
    another singular data (e.g., a Nino index)."""
    data = eval(f"gridded_da_{time_type}")()
    x = data.isel(lat=0, lon=0)
    y = data if gridded else data.isel(lat=0, lon=0)
    if single_arg:
        result = func(y, dim="time")
    else:
        result = func(x, y, "time")
    assert not result.isnull().any()
    if func == linear_slope:
        expected_shape = (3, 3) if gridded else ()
    elif func == linregress:
        expected_shape = (3, 3, 5) if gridded else (5,)
    assert result.shape == expected_shape

@pytest.mark.parametrize("time_type", TIME_TYPES)
@pytest.mark.parametrize("func", (linear_slope, linregress))
def test_linear_regression_gridded_to_gridded(
    gridded_da_datetime, gridded_da_cftime, gridded_da_float, time_type, func
):
    """Tests that ``linear_slope`` works between one grid of data and another grid
    of data."""
    data1 = eval(f"gridded_da_{time_type}")()
    data2 = eval(f"gridded_da_{time_type}")()
    result = func(data1, data2, "time")
    assert not result.isnull().any()
    expected_shape = (3, 3) if func == linear_slope else (3, 3, 5)
    assert result.shape == expected_shape

@pytest.mark.parametrize("gridded", (True, False))
@pytest.mark.parametrize("time_type", TIME_TYPES)
@pytest.mark.parametrize("func", (linear_slope, linregress))
@pytest.mark.parametrize("single_arg", (True, False))
def test_linear_regression_time_ds(
    gridded_ds_datetime,
    gridded_ds_cftime,
    gridded_ds_float,
    time_type,
    gridded,
    func,
    single_arg,
):
    """Tests that linear slope can be computed on an in-memory Dataset with respect to
    time."""
    ds = eval(f"gridded_ds_{time_type}")
    x = ds["time"]
    if not gridded:
        ds = ds.isel(lat=0, lon=0)
    if single_arg:
        result = func(ds, dim="time")
    else:
        result = func(x, ds, "time")
    assert len(ds.data_vars) > 1
    assert not result["foo"].isnull().any()
    assert not result["bar"].isnull().any()
    if func == linear_slope:
        expected_shape = (3, 3) if gridded else ()
    elif func == linregress:
        expected_shape = (3, 3, 5) if gridded else (5,)
    assert result["foo"].shape == expected_shape
    assert result["bar"].shape == expected_shape

@pytest.mark.parametrize("gridded", (True, False))
@pytest.mark.parametrize("time_type", TIME_TYPES)
@pytest.mark.parametrize("func", (linear_slope, linregress))
@pytest.mark.parametrize("single_arg", (True, False))
def test_linear_regression_da_dask(
    gridded_da_datetime,
    gridded_da_cftime,
    gridded_da_float,
    time_type,
    gridded,
    func,
    single_arg,
):
    """Tests that linear slope can be computed on a dask DataArray."""
    data = eval(f"gridded_da_{time_type}")().chunk()
    x = data["time"]
    y = data if gridded else data.isel(lat=0, lon=0)
    if single_arg:
        expected = func(y.load(), dim="time")
        actual = func(y, dim="time").compute()
    else:
        expected = func(x.load(), y.load(), "time")
        actual = func(x, y, "time").compute()
    assert_allclose(expected, actual)

@pytest.mark.parametrize("gridded", (True, False))
@pytest.mark.parametrize("time_type", TIME_TYPES)
@pytest.mark.parametrize("func", (linear_slope, linregress))
@pytest.mark.parametrize("single_arg", (True, False))
def test_linear_slope_ds_dask(
    gridded_ds_datetime,
    gridded_ds_cftime,
    gridded_ds_float,
    time_type,
    gridded,
    func,
    single_arg,
):
    """Tests that linear slope can be computed on a dask Dataset."""
    ds = eval(f"gridded_ds_{time_type}").chunk()
    x = ds["time"]
    if not gridded:
        ds = ds.isel(lat=0, lon=0)
    if single_arg:
        expected = func(ds.load(), dim="time")
        actual = func(ds, dim="time").compute()
    else:
        expected = func(x.load(), ds.load(), "time")
        actual = func(x, ds, "time").compute()
    assert_allclose(expected, actual)

def test_linear_slope_time_conversion_accurate(ts_annual_da):
    """Tests that computed slope with cftime/datetime conversion is accurate."""
    data = ts_annual_da()
    data2 = data.copy()
    # Just create numeric time, e.g. time in years.
    data2["time"] = np.arange(data2.time.size)
    slope_cftime = linear_slope(data["time"], data, "time")
    slope_ints = linear_slope(data2["time"], data2, "time")
    diff = slope_cftime - slope_ints
    assert np.abs(diff) < 1e-4

def test_linear_slope_raises_error(ts_annual_da):
    """Tests that error is raised for `linear_slope` when independent variable is
    in y position."""
    with pytest.raises(ValueError):
        data = ts_annual_da()
        x = data["time"]
        y = data
        linear_slope(y, x, "time")

def test_linear_slope_linregress_same_result(ts_annual_da):
    """Tests that ``linear_slope`` and ``linregress`` returns the same result
    for the slope parameter."""
    data = ts_annual_da()
    x = data["time"]
    y = data
    m1 = linear_slope(x, y, "time")
    m2 = linregress(x, y, "time").sel(parameter="slope")
    assert np.abs(m1.values - m2.values) < 1e-15

max_stars_repo_name: luke-gregor/esmtools

max_issues_repo_name: luke-gregor/esmtools

max_forks_repo_name: luke-gregor/esmtools


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

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

import tensorflow as tf
from sagemakercv.core import training_ops, GenericRoIExtractor, TargetEncoder, BoxDetector, RandomSampler
from ..builder import HEADS, build_box_head, build_mask_head

class StandardRoIHead(tf.keras.Model):
    """Simplest base roi head including one bbox head and one mask head."""

    def __init__(self,
                 bbox_head,
                 bbox_roi_extractor,
                 bbox_sampler,
                 box_encoder,
                 inference_detector,
                 mask_head=None,
                 mask_roi_extractor=None,
                 name="StandardRoIHead",
                 trainable=True,
                 *args,
                 **kwargs,
                 ):
        super(StandardRoIHead, self).__init__(name=name, trainable=trainable, *args, **kwargs)
        self.bbox_head = bbox_head
        self.bbox_roi_extractor = bbox_roi_extractor
        self.bbox_sampler = bbox_sampler
        self.box_encoder = box_encoder
        self.mask_head = mask_head
        self.mask_roi_extractor = mask_roi_extractor if mask_roi_extractor is not None \
                                  else self.bbox_roi_extractor
        self.inference_detector = inference_detector

    @property
    def with_mask(self):
        """bool: whether the RoI head contains a `mask_head`"""
        return hasattr(self, 'mask_head') and self.mask_head is not None

    def call(self,
             fpn_feats,
             img_info,
             proposals,
             gt_bboxes=None,
             gt_labels=None,
             gt_masks=None,
             training=True):
        model_outputs=dict()

        if training:
            box_targets, class_targets, rpn_box_rois, proposal_to_label_map = self.bbox_sampler(proposals,
                                                                                gt_bboxes,
                                                                                gt_labels)
        else:
            rpn_box_rois = proposals

        box_roi_features = self.bbox_roi_extractor(fpn_feats, rpn_box_rois)

        class_outputs, box_outputs, _ = self.bbox_head(inputs=box_roi_features)

        if not training:
            model_outputs.update(self.inference_detector(class_outputs, box_outputs, rpn_box_rois, img_info))
            model_outputs.update({'class_outputs': tf.nn.softmax(class_outputs),
                                  'box_outputs': box_outputs,
                                  'anchor_boxes': rpn_box_rois})
        else:
            if self.bbox_head.loss.box_loss_type not in ["giou", "ciou"]:
                encoded_box_targets = self.box_encoder(boxes=rpn_box_rois,
                                                       gt_bboxes=box_targets,
                                                       gt_labels=class_targets)
            model_outputs.update({
                'class_outputs': class_outputs,
                'box_outputs': box_outputs,
                'class_targets': class_targets,
                'box_targets': encoded_box_targets if self.bbox_head.loss.box_loss_type \
                                                      not in ["giou", "ciou"] \
                                                      else box_targets,
                'box_rois': rpn_box_rois,
            })
            total_loss, class_loss, box_loss = self.bbox_head.loss(model_outputs['class_outputs'],
                                                                   model_outputs['box_outputs'],
                                                                   model_outputs['class_targets'],
                                                                   model_outputs['box_targets'],
                                                                   model_outputs['box_rois'],
                                                                   img_info)
            model_outputs.update({
                'total_loss_bbox': total_loss,
                'class_loss': class_loss,
                'box_loss': box_loss
            })
        if not self.with_mask:
            return model_outputs
        if not training:
            return self.call_mask(model_outputs, fpn_feats, training=False)
        max_fg = int(self.bbox_sampler.batch_size_per_im * self.bbox_sampler.fg_fraction)
        return self.call_mask(model_outputs,
                              fpn_feats,
                              class_targets=class_targets,
                              box_targets=box_targets,
                              rpn_box_rois=rpn_box_rois,
                              proposal_to_label_map=proposal_to_label_map,
                              gt_masks=gt_masks,
                              max_fg=max_fg,
                              training=True)

    def call_mask(self,
                  model_outputs,
                  fpn_feats,
                  class_targets=None,
                  box_targets=None,
                  rpn_box_rois=None,
                  proposal_to_label_map=None,
                  gt_masks=None,
                  max_fg=None,
                  training=True):
        if not training:
            selected_box_rois = model_outputs['detection_boxes']
            class_indices = model_outputs['detection_classes']
            class_indices = tf.cast(class_indices, dtype=tf.int32)

        else:
            selected_class_targets, selected_box_targets, \
            selected_box_rois, proposal_to_label_map = training_ops.select_fg_for_masks(
                class_targets=class_targets,
                box_targets=box_targets,
                boxes=rpn_box_rois,
                proposal_to_label_map=proposal_to_label_map,
                max_num_fg=max_fg
            )

            class_indices = selected_class_targets
            class_indices = tf.cast(selected_class_targets, dtype=tf.int32)

        mask_roi_features = self.mask_roi_extractor(
                fpn_feats,
                selected_box_rois,
            )

        mask_outputs = self.mask_head(inputs=mask_roi_features, class_indices=class_indices)

        if training:
            mask_targets = training_ops.get_mask_targets(
                fg_boxes=selected_box_rois,
                fg_proposal_to_label_map=proposal_to_label_map,
                fg_box_targets=selected_box_targets,
                mask_gt_labels=gt_masks,
                output_size=self.mask_head._mrcnn_resolution
            )

            model_outputs.update({
                'mask_outputs': mask_outputs,
                'mask_targets': mask_targets,
                'selected_class_targets': selected_class_targets,
            })

            mask_loss = self.mask_head.loss(model_outputs['mask_outputs'],
                                             model_outputs['mask_targets'],
                                             model_outputs['selected_class_targets'],)
            model_outputs.update({'mask_loss': mask_loss})

        else:
            model_outputs.update({
                'detection_masks': tf.nn.sigmoid(mask_outputs),
            })

        return model_outputs

@HEADS.register("StandardRoIHead")
def build_standard_roi_head(cfg):
    roi_head = StandardRoIHead
    bbox_head = build_box_head(cfg)
    bbox_roi_extractor = GenericRoIExtractor(cfg.MODEL.FRCNN.ROI_SIZE,
                                            cfg.MODEL.FRCNN.GPU_INFERENCE)
    bbox_sampler = RandomSampler(batch_size_per_im=cfg.MODEL.RCNN.BATCH_SIZE_PER_IMAGE,
                                 fg_fraction=cfg.MODEL.RCNN.FG_FRACTION,
                                 fg_thresh=cfg.MODEL.RCNN.THRESH,
                                 bg_thresh_hi=cfg.MODEL.RCNN.THRESH_HI,
                                 bg_thresh_lo=cfg.MODEL.RCNN.THRESH_LO)
    box_encoder = TargetEncoder(bbox_reg_weights=cfg.MODEL.BBOX_REG_WEIGHTS)
    inference_detector = BoxDetector(use_batched_nms=cfg.MODEL.INFERENCE.USE_BATCHED_NMS,
                                     rpn_post_nms_topn=cfg.MODEL.INFERENCE.POST_NMS_TOPN,
                                     detections_per_image=cfg.MODEL.INFERENCE.DETECTIONS_PER_IMAGE,
                                     test_nms=cfg.MODEL.INFERENCE.DETECTOR_NMS,
                                     class_agnostic_box=cfg.MODEL.INFERENCE.CLASS_AGNOSTIC,
                                     bbox_reg_weights=cfg.MODEL.BBOX_REG_WEIGHTS)

    if cfg.MODEL.INCLUDE_MASK:
        mask_roi_extractor = GenericRoIExtractor(cfg.MODEL.MRCNN.ROI_SIZE,
                                                 cfg.MODEL.MRCNN.GPU_INFERENCE)
        mask_head = build_mask_head(cfg)
    else:
        mask_head = None
        mask_roi_extractor = None
    return roi_head(bbox_head=bbox_head,
                    bbox_roi_extractor=bbox_roi_extractor,
                    bbox_sampler=bbox_sampler,
                    box_encoder=box_encoder,
                    inference_detector=inference_detector,
                    mask_head=mask_head,
                    mask_roi_extractor=mask_roi_extractor)

max_stars_repo_name: aws-samples/amazon-sagemaker-cv

max_issues_repo_name: aws-samples/amazon-sagemaker-cv

max_forks_repo_name: aws-samples/amazon-sagemaker-cv

reason: from sage


text:
import numpy as np
import matplotlib.pyplot as plt
from random import randint, random
from numpy.core.fromnumeric import diagonal

'''
Load data from file
'''
def load_data(filename):
    file = open(filename, 'r')
    tmp_str = file.readline()
    tmp_arr = tmp_str[:-1].split(' ')
    N = int(tmp_arr[0])
    n_row = int(tmp_arr[1])
    n_col = int(tmp_arr[2])
    print('N=%d, row=%d, col=%d' %(N,n_row,n_col))
    data = np.zeros([N, n_row * n_col + 1])
    for n in range(N):
        tmp_str = file.readline()
        tmp_arr = tmp_str[:-1].split(' ')
        for i in range(n_row * n_col + 1):
            data[n][i] = int(tmp_arr[i])
    file.close()
    return N, n_row, n_col, data

'''
Sigmoid function
'''
def sigmoid(s):
    large=30
    if s<-large: s=-large
    if s>large: s=large
    return (1 / (1 + np.exp(-s)))

'''
Cost/Loss funtion
'''
def cost(X, Y, N, W, b, v, c):
    epsi = 1.e-12
    sum = 0
    for n in range(N):
        prediction = predict(X[n], J, W, b, v, c)
        # verify if the value of the prediction is in [1e-12, 1-1e-12]
        if prediction < epsi:
            prediction = epsi
        if prediction > 1 - epsi:
            prediction = 1 - epsi
        sum += Y[n] * np.log(prediction) + (1 - Y[n]) * np.log(1 - prediction)
    E = - sum / N
    return E

'''
Calculate h
'''
def get_h(x, W, b):
    # x, w, b -> b + w * x -> s
    s = b + np.dot(W, x)
    # s -> s.map(sigmoid) -> h
    func = lambda a: sigmoid(a)
    h = np.array([func(s_j) for s_j in s])
    #print('a')
    return h

'''
Predict label
'''
def predict(x, J, W, b, v, c):
    h = get_h(x, W, b)
    # c, v, h -> c + v.T * h -> z
    sum = 0
    for j in range(J):
        sum += v[j] * h[j]
    z = c + sum
    # prediction
    y = sigmoid(z)
    return y

'''
Calculate the gradients for w, b, v, c
'''
def gradient(x, y, v, h, prediction):

    func = lambda a: a * (1 - a)
    diagonal = np.array([func(h_m) for h_m in h])
    matrix1 = np.diag(diagonal) # J x J
    matrix2 = np.tile(x, (J, 1))
    matrix2 = np.transpose(matrix2) * v
    matrix2 = np.transpose(matrix2) # J x I
    delta = prediction - y
    Gw = delta * np.matmul(matrix1, matrix2) # w's gradient
    Gb = delta * np.dot(matrix1, v)          # b's gradient
    Gv = delta * h                           # v's gradient
    Gc = delta                               # c's gradient
    '''
    delta = prediction - y
    diagonal = np.multiply(h, (1 - h))
    Gb = delta * np.multiply(v, diagonal)
    Gw = np.tensordot(Gb, x, axes=0)
    Gv = delta * h
    Gc = delta
    '''
    return Gw, Gb, Gv, Gc

'''
Update w, b, v, c
'''
def update(x, y, J, eta, W, b, v, c):
    # calculate the gradients
    h = get_h(x, W, b)
    prediction = predict(x, J, W, b, v, c)
    Gw, Gb, Gv, Gc = gradient(x, y, v, h, prediction)
    # update the first layer
    W = W - (eta * Gw) # w (t + 1)
    b = b - (eta * Gb) # b (t + 1)
    # update the second layer
    v = v - (eta * Gv) # v (t + 1)
    c = c - (eta * Gc) # c (t + 1)
    return W, b, v, c

'''
Run shallow logistic classifier
'''
def run_slc(X, Y, N, J, subloop, eta, max_iteration, w, b, v, c, errors):
    #epsi = 10e-3
    epsi = 0
    iteration = 0
    while (errors[-1] > epsi):
        for j in range(subloop):
            # choose random data from dataset
            m = randint(0, N - 1)
            x = X[m]
            y = Y[m]
            # update w, b, v, c from (t) to (t + 1)
            w, b, v, c = update(x, y, J, eta, w, b, v, c)
            # calculate error
        error = cost(X, Y, N, w, b, v, c)
        errors.append(error)
        print('number of iterations = %d, cost = %f, eta = %e\r' %(iteration, errors[-1], eta), end='')
        iteration += subloop
        if (iteration > max_iteration):
            break
    return w, b, v, c, errors

def plot_error(err):
    plt.plot(range(len(err)), err, marker='o')
    plt.xlabel('Iterations')
    plt.ylabel('Number of misclassifications')
    plt.ylim([0,5])
    plt.show()
    return

def confusion(Xeval,Yeval,N, J, w, b, v, c):
    C=np.zeros([2,2])
    for n in range(N):
        y = predict(Xeval[n], J, w, b, v, c)
        if(y<0.5 and Yeval[n]<0.5): C[0,0]=C[0,0]+1
        if(y>0.5 and Yeval[n]>0.5): C[1,1]=C[1,1]+1
        if(y<0.5 and Yeval[n]>0.5): C[1,0]=C[1,0]+1
        if(y>0.5 and Yeval[n]<0.5): C[0,1]=C[0,1]+1
    return C

'''
Main execution
'''

# load data from dataset
#dataset = './Data/XOR.txt'
#dataset = './Data/line600.txt'      # almost linearly separable (J = 2)
dataset = './Data/rectangle600.txt' # it's harder to separate linearly (J = 10)
#dataset = './Data/square_circle.txt' # linearly separable (J = 1)
N, n_row, n_col, data = load_data(dataset)
I = n_row * n_col

# shuffle dataset
np.random.shuffle(data)

# split data into training and testing sets
N_train = int(N * 0.8)
X_train = data[:N_train, :-1]
Y_train = data[:N_train, -1]
N_test = N - N_train
X_test = data[N_train:N, :-1]
Y_test = data[N_train:N, -1]

# replace labels -1 for 0
np.place(Y_train, Y_train != 1, [0])
np.place(Y_test, Y_test != 1, [0])

J = 10 # number of neurons in the hidden layer

# initialize w, b, v and c with random values
W = np.random.randn(J, I)
b = np.random.randn(J)
v = np.random.rand(J)
c = np.random.rand(1)

eta = 1 # learning rate

Nbiter = 3000 # number of iterations
subloop = 300 # only calculate the cost in each 100 iterations

# initialize errors' list
errors = []
errors.append(cost(X_train, Y_train, N_train, W, b, v, c))

W, b, v, c, errors = run_slc(X_train, Y_train, N_train, J, subloop, eta, Nbiter, W, b, v, c, errors); print("\n"); eta = 0.5 * eta
W, b, v, c, errors = run_slc(X_train, Y_train, N_train, J, subloop, eta, Nbiter, W, b, v, c, errors); print("\n"); eta = 0.5 * eta
#W, b, v, c, errors = run_slc(X_train, Y_train, N_train, J, subloop, eta, Nbiter, W, b, v, c, errors); print("\n"); eta = 0.5 * eta
#W, b, v, c, errors = run_slc(X_train, Y_train, N_train, J, subloop, eta, Nbiter, W, b, v, c, errors); print("\n")

plot_error(errors)

print('in-samples error = %f ' % (cost(X_train, Y_train, N_train, W, b, v, c)))
# calculate and show confusion matrix
C = confusion(X_train, Y_train, N_train, J, W, b, v, c)
print(C)

print('out-samples error=%f' % (cost(X_test, Y_test, N_test, W, b, v, c)))
# calculate and show confusion matrix
C =confusion(X_test, Y_test, N_test, J, W, b, v, c)
print(C)

max_stars_repo_name: PedromfRibeiro/OML-SLC

max_issues_repo_name: PedromfRibeiro/OML-SLC

max_forks_repo_name: PedromfRibeiro/OML-SLC

reason: import numpy,from numpy


text:
"""
Tests whether the fluent API works
"""
import os

import networkx
from testfixtures import compare

from mlinspect.inspections.materialize_first_rows_inspection import MaterializeFirstRowsInspection
from mlinspect.utils import get_project_root
from mlinspect.pipeline_inspector import PipelineInspector
from .utils import get_expected_dag_adult_easy_ipynb, get_expected_dag_adult_easy_py

ADULT_EASY_FILE_PY = os.path.join(str(get_project_root()), "test", "pipelines", "adult_easy.py")
FILE_NB = os.path.join(str(get_project_root()), "test", "pipelines", "adult_easy.ipynb")

def test_inspector_adult_easy_py_pipeline():
    """
    Tests whether the .py version of the inspector works
    """
    inspection_result = PipelineInspector\
        .on_pipeline_from_py_file(ADULT_EASY_FILE_PY)\
        .add_inspection(MaterializeFirstRowsInspection(5))\
        .execute()
    extracted_dag = inspection_result.dag
    expected_dag = get_expected_dag_adult_easy_py()
    compare(networkx.to_dict_of_dicts(extracted_dag), networkx.to_dict_of_dicts(expected_dag))

def test_inspector_adult_easy_ipynb_pipeline():
    """
    Tests whether the .ipynb version of the inspector works
    """
    inspection_result = PipelineInspector\
        .on_pipeline_from_ipynb_file(FILE_NB)\
        .add_inspection(MaterializeFirstRowsInspection(5))\
        .execute()
    extracted_dag = inspection_result.dag
    expected_dag = get_expected_dag_adult_easy_ipynb()
    compare(networkx.to_dict_of_dicts(extracted_dag), networkx.to_dict_of_dicts(expected_dag))

def test_inspector_adult_easy_str_pipeline():
    """
    Tests whether the str version of the inspector works
    """
    with open(ADULT_EASY_FILE_PY) as file:
        code = file.read()

        inspection_result = PipelineInspector\
            .on_pipeline_from_string(code)\
            .add_inspection(MaterializeFirstRowsInspection(5))\
            .execute()
        extracted_dag = inspection_result.dag
        expected_dag = get_expected_dag_adult_easy_py()
        assert networkx.to_dict_of_dicts(extracted_dag) == networkx.to_dict_of_dicts(expected_dag)

max_stars_repo_name: stefan-grafberger/mlinspect-cidr

max_issues_repo_name: stefan-grafberger/mlinspect-cidr

max_forks_repo_name: stefan-grafberger/mlinspect-cidr

reason: import networkx


text:
# coding=utf-8
#通过指定索引标签读取数据
import numpy as np
np.random.seed(2018)
nd11=np.random.random([10])
print(nd11)
print(nd11[3])          #获取第3+1个元素
print(nd11[3: 6])        #获取第3+1到第6个元素
print(nd11[1: 6: 2])      #从第1+1到第6，步长为2获取元素
print(nd11[: : 2])        #从第1开始，步长为2获取元素
nd12=np.arange(25).reshape([5,5])           #此处reshape将原矩阵按顺序重塑为5行5列矩阵
print((nd12))
print(nd12[1: 3, 1: 3])                        #截取2、3行，2、3列元素
print(nd12[(nd12 > 3) & (nd12 < 10)])
print(nd12[[1, 2]])                          #此处截取2、3行，原理不明，容易理解的方式是print(nd12[1:3,:])
print(nd12[:, 1: 3])

'''输出
[0.88234931 0.10432774 0.90700933 0.3063989  0.44640887 0.58998539
 0.8371111  0.69780061 0.80280284 0.10721508]
0.3063988986063515
[0.3063989  0.44640887 0.58998539]
[0.10432774 0.3063989  0.58998539]
[0.88234931 0.90700933 0.44640887 0.8371111  0.80280284]
[[ 0  1  2  3  4]
 [ 5  6  7  8  9]
 [10 11 12 13 14]
 [15 16 17 18 19]
 [20 21 22 23 24]]
[[ 6  7]
 [11 12]]
[4 5 6 7 8 9]
[[ 5  6  7  8  9]
 [10 11 12 13 14]]
[[ 1  2]
 [ 6  7]
 [11 12]
 [16 17]
 [21 22]]
'''


text:
import numpy as np
from hypothesis import given, settings, strategies as st

import est_dir

@settings(max_examples=5, deadline=None)
@given(st.integers(2, 100), st.integers(1, 10), st.integers(16, 100))
def test_1(m, min_eig, max_eig):
    """
    Series of tests for quad_func_params() to check oututs are of correct
    form when max_eig != min_eig.
    """
    mat = est_dir.quad_func_params(min_eig, max_eig, m)
    assert(mat.shape == (m, m))
    for j in range(m):
        for i in range(m):
            if i != j:
                assert(np.round(mat[j, i], 4) == np.round(mat[j, i], 4))

@settings(max_examples=5, deadline=None)
@given(st.integers(2, 500))
def test_2(m, ):
    """
    Series of tests for quad_func_params() to check oututs are of correct
    form when max_eig == min_eig.
    """
    min_eig = 1
    max_eig = 1
    A = est_dir.quad_func_params(min_eig, max_eig, m)
    assert(np.all(np.round(A, 3) == np.identity(m)))

max_stars_repo_name: Megscammell/Estimating-Direction

max_issues_repo_name: Megscammell/Estimating-Direction

max_forks_repo_name: Megscammell/Estimating-Direction


text:
# -*- coding: utf-8 -*-
"""
@author: Dung D. Le (Andrew) <ddle.2015@smu.edu.sg>
"""

import numpy as np
import random
from ...utils.data_utils import Dataset

try:
    import torch
except ImportError:
    torch = None

"""Firstly, we define a helper function to generate\sample training ordinal triplets:
   Step 1:
   given rated item i, randomly choose item j and check whether rating of j is lower than i,
   if not randomly sample another item.
   each row of the sampled data in the following form:
        [userId itemId_i itemId_j]
   for each user u, he/she prefers item i over item j.
   """

def sample_data(X, data):
    sampled_data = np.zeros((data.shape[0], 5), dtype=np.int)
    data = data.astype(int)

    for k in range(0, data.shape[0]):
        u = data[k, 0]
        i = data[k, 1]
        ratingi = data[k, 2]
        j = random.randint(0, X.shape[1] - 1)

        while X[u, j] > ratingi:
            j = random.randint(0, X.shape[1] - 1)

        sampled_data[k, :] = [u, i, j, ratingi, X[u, j]]

    return sampled_data

def ibpr(X, data, k, lamda=0.001, n_epochs=150, learning_rate=0.05, batch_size=100, init_params=None):
    # X = sp.csr_matrix(X)
    Data = Dataset(data)

    # Initial user factors
    if init_params['U'] is None:
        U = torch.randn(X.shape[0], k, requires_grad=True)
    else:
        U = init_params['U']
        U = torch.from_numpy(U)

    # Initial item factors
    if init_params['V'] is None:
        V = torch.randn(X.shape[1], k, requires_grad=True)
    else:
        V = init_params['V']
        V = torch.from_numpy(V)

    optimizer = torch.optim.Adam([U, V], lr=learning_rate)
    for epoch in range(n_epochs):
        num_steps = int(Data.data.shape[0] / batch_size)
        for i in range(1, num_steps + 1):
            batch_c, _ = Data.next_batch(batch_size)
            # print(batch_c, idx)
            sampled_batch = sample_data(X, batch_c)

            regU = U[sampled_batch[:, 0], :]
            regI = V[sampled_batch[:, 1], :]
            regJ = V[sampled_batch[:, 2], :]

            regU_unq = U[np.unique(sampled_batch[:, 0]), :]
            regI_unq = V[np.unique(sampled_batch[:, 1:]), :]

            regU_norm = regU / regU.norm(dim=1)[:, None]
            regI_norm = regI / regI.norm(dim=1)[:, None]
            regJ_norm = regJ / regJ.norm(dim=1)[:, None]

            Scorei = torch.acos(torch.clamp(torch.sum(regU_norm * regI_norm, dim=1), -1 + 1e-7, 1 - 1e-7))
            Scorej = torch.acos(torch.clamp(torch.sum(regU_norm * regJ_norm, dim=1), -1 + 1e-7, 1 - 1e-7))

            loss = lamda * (regU_unq.norm().pow(2) + regI_unq.norm().pow(2)) - torch.log(
                torch.sigmoid(Scorej - Scorei)).sum()
            optimizer.zero_grad()
            loss.backward()
            optimizer.step()
        print('epoch:', epoch, 'loss:', loss)

    # since the user's preference is defined by the angular distance, we can normalize the user/item vectors without changing the ranking
    U = torch.nn.functional.normalize(U, p=2, dim=1)
    V = torch.nn.functional.normalize(V, p=2, dim=1)
    U = U.data.cpu().numpy()
    V = V.data.cpu().numpy()

    res = {'U': U, 'V': V}

    return res

max_stars_repo_name: GuoJingyao/cornac

max_issues_repo_name: GuoJingyao/cornac

max_forks_repo_name: GuoJingyao/cornac


text:
import cv2 as cv
import numpy as np
from matplotlib import pyplot as plt

img = cv.imread('gradient.jpg', 0)
_, th1 = cv.threshold(img, 50, 255, cv.THRESH_BINARY)
_, th2 = cv.threshold(img, 200, 255, cv.THRESH_BINARY_INV)
_, th3 = cv.threshold(img, 127, 255, cv.THRESH_TRUNC)
_, th4 = cv.threshold(img, 127, 255, cv.THRESH_TOZERO)
_, th5 = cv.threshold(img, 127, 255, cv.THRESH_TOZERO_INV)

titles = ['Original Image', 'BINARY', 'BINARY_INV', 'TRUNC', 'TOZERO', 'TOZERO_INV']
images = [img, th1, th2, th3, th4, th5]

for i in range(6):
    #  show as (2,3) dimension
    plt.subplot(2, 3, i+1), plt.imshow(images[i], 'gray')
    #  title in each image
    plt.title(titles[i])
    #  to hide (x, y) dimension graph
    plt.xticks([]), plt.yticks([])

plt.show()

max_stars_repo_name: KyawThuHtun/OpenCV-with-Python

max_issues_repo_name: KyawThuHtun/OpenCV-with-Python

max_forks_repo_name: KyawThuHtun/OpenCV-with-Python


text:
#!/usr/bin/env python
'''
  summarise numerical columns
'''

import argparse
import collections
import csv
import logging
import sys

import numpy

def check(summary, col):
  if col not in summary:
    summary[col] = {'n': 0, 'sum': 0, 'max': sys.float_info.min, 'min': sys.float_info.max, 'd': []}
  return summary

def main(delimiter, x, y, z, fh, out):
  logging.info('reading...')
  cell = {}
  xs = set()
  ys = set()
  for row in csv.DictReader(fh, delimiter=delimiter):
    xs.add(row[x])
    ys.add(row[y])
    cell[(row[x], row[y])] = row[z]

  fo = csv.DictWriter(out, delimiter=delimiter, fieldnames=[x] + sorted(list(xs)))
  fo.writeheader()
  for y in sorted(list(ys)):
    row = {x: y}
    for xv in xs:
      row[xv] = cell.get((xv, y),'')
    fo.writerow(row)

if __name__ == '__main__':
  parser = argparse.ArgumentParser(description='Basic stats of specified columns')
  parser.add_argument('--delimiter', required=False, default=',', help='input files')
  parser.add_argument('--x', required=True, help='x column')
  parser.add_argument('--y', required=True, help='y column')
  parser.add_argument('--z', required=True, help='z column')
  parser.add_argument('--verbose', action='store_true', help='more logging')
  parser.add_argument('--encoding', default='utf-8', help='file encoding')
  parser.add_argument('--quiet', action='store_true', help='more logging')
  args = parser.parse_args()
  if args.verbose:
    logging.basicConfig(format='%(asctime)s %(levelname)s %(message)s', level=logging.DEBUG)
  elif args.quiet:
    logging.basicConfig(format='%(asctime)s %(levelname)s %(message)s', level=logging.WARN)
  else:
    logging.basicConfig(format='%(asctime)s %(levelname)s %(message)s', level=logging.INFO)

  if "reconfigure" in dir(sys.stdin):
    sys.stdin.reconfigure(encoding=args.encoding)
  main(args.delimiter, args.x, args.y, args.z, sys.stdin, sys.stdout)

max_stars_repo_name: supernifty/csvtools

max_issues_repo_name: supernifty/csvtools

max_forks_repo_name: supernifty/csvtools


text:
import numpy as np

def djisktra_slow():
 #Initialize largest value of shortest distance to large number
 maxDist = 1000000
 #Load graph from file and populate Path weight matrix
 #Note that if weight for a path is equal to maxDist it means given path is not connected
 f = open("djikstra.txt","r")
 lines = [line.rstrip('\n').split() for line in f]
 Lij = maxDist*np.ones((len(lines),len(lines)),dtype=int)
 for g in lines:
    for j in g:
        list = j.split(',')
        if len(list) > 1:
           Lij[lines.index(g)][int(list[0])-1] = int(list[1])
 f.close()

 #Initialize full graph (use set so that comparison with X can be done irrespective of order
 #in which nodes are captured by algorithm)
 V = {i for i in range(1,len(lines)+1)}
 # Initialize portion of graph which Djikstra algorithm has already seen
 X = {1}
 A = {}
 A[1] = 0
 # Initialize counters
 iteration = 0
 cnt = 0

 while X != V:
    LowDist = maxDist
    LowIdx = 0
    for i in X:
        for j in V:
            if j not in X:
               dist = A[i] + Lij[i-1][j-1]
               cnt += 1
               #print(i,j,dist)
               if dist < LowDist:
                   LowDist = dist
                   LowIdx = j

    X.add(LowIdx)
    A[LowIdx] = LowDist
    iteration += 1

 print(iteration,cnt)
 print(A[7],A[37],A[59],A[82],A[99],A[115],A[133],A[165],A[188],A[197])

djisktra_slow()

#1 - There is a shortest path s-t with no repeated vertices; Shortest path may have as many as n-1 vertices
#2 - When all edges are distinct powers of 2
#3 - Always
#4 - Always terminates but in some cases paths computed will not be shortest; always terminates and in some cases paths computed will be correct shortest path
#5 - Djkistra algorithm always terminates


text:
#!/usr/bin/env python3

import cheminfo

import sys

import numpy as np
import qml
import scipy
from qml.kernels import (get_atomic_local_gradient_kernel,
                         get_atomic_local_kernel)
from qml.math import svd_solve
from qml.representations import generate_fchl_acsf
import rmsd
np.random.seed(666)

FILENAME_TEST = "data/ethanol_ccsd_t-test.npz"
FILENAME_TRAIN = "data/ethanol_ccsd_t-train.npz"

def get_data_from_file(filename, n=100):

    data = np.load(filename)

    X  = []
    dX = []
    Q  = []

    E  = []
    F  = []

    max_n = len(data["E"])

    index = np.random.choice(max_n, size=n, replace=True)

    nuclear_charges = data["z"]
    max_atoms = len(nuclear_charges)

    for i in index:

        coordinates = data["R"][i]

        (rep, drep) = generate_fchl_acsf(nuclear_charges, coordinates,
               gradients=True, pad=max_atoms)

        X.append(rep)
        dX.append(drep)
        Q.append(nuclear_charges)
        E.append(data["E"][i])
        F.append(data["F"][i])

    X  = np.array(X)
    dX = np.array(dX)
    E  = np.array(E).flatten()
    F  = np.array(F)

    return X, dX, Q, E, F

def train(dataname, n_train=100):

    SIGMA = 10.0

    filename_train = "data/" + dataname + "-train.npz"

    # Read training data from file
    X,  dX,  Q,  E,  F  = get_data_from_file(filename_train, n=n_train)

    offset = E.mean()
    E  -= offset
    print("OFFSET: ", offset)

    F = np.concatenate(F)
    Y = np.concatenate((E, F.flatten()))

    print("Generating Kernels ...")
    Kte = get_atomic_local_kernel(X, X,  Q, Q,  SIGMA)
    Kt = get_atomic_local_gradient_kernel(X, X,  dX,  Q, Q,  SIGMA)

    C = np.concatenate((Kte, Kt))

    print("Alphas operator ...")
    alpha = svd_solve(C, Y, rcond=1e-11)

    np.save("data/"+dataname+"_offset.npy", offset)
    np.save("data/"+dataname+"_sigma.npy", SIGMA)
    np.save("data/"+dataname+"_alphas.npy", alpha)
    np.save("data/"+dataname+"_Q.npy", Q)
    np.save("data/"+dataname+"_X.npy", X)

    return

if __name__ == "__main__":

    import argparse
    import sys

    parser = argparse.ArgumentParser()
    parser.add_argument('--correct', action='store_true', help='Add penalty and correction molecules to kernel')
    parser.add_argument('dataname', metavar='NAME', type=str)

    args = parser.parse_args()

    if len(sys.argv[1:]) < 1:
        print("choose some data")
        quit()

    dataname = args.dataname

    train(dataname, n_train=150)

    # Get sample from test
    data = np.load("data/" + dataname + "-test.npz")
    coordinates = data["R"][0]
    nuclear_charges = data["z"]
    nuclear_charges = [cheminfo.convert(x) for x in nuclear_charges]

    txt = rmsd.set_coordinates(nuclear_charges, coordinates)
    f = open("examples/"+dataname+".xyz", 'w')
    f.write(txt)
    f.close()

max_stars_repo_name: andersx/qml-ase

max_issues_repo_name: andersx/qml-ase

max_forks_repo_name: andersx/qml-ase

reason: import numpy,import scipy


text:
# Licensed with the 3-clause BSD license.  See LICENSE for details.
import pytest
from logging import Logger
from itertools import repeat

import numpy as np
import astropy.units as u
import sqlalchemy as sa
from sbpy.data import Orbit

from .. import util, schema
from ..schema import Obj, Obs, Eph, Found
from ..util import RADec, FieldOfView
from ..db import SBDB
from ..exceptions import BadObjectID, NoEphemerisError, SourceNotFoundError
from .skytiles import N_tiles

@pytest.fixture
def db():
    db = SBDB.create_test_db()
    yield db
    db.close()

class TestSBDB:
    def test_verify_database(self, db):
        logger = Logger('test')
        db.verify_database(logger)

        metadata = sa.MetaData(bind=db.engine)
        metadata.reflect()
        for table in db.DB_NAMES:
            assert table in metadata.tables.keys()

        db.session.execute('DROP TABLE eph')
        metadata.clear()
        metadata.reflect()
        assert 'eph' not in metadata.tables.keys()
        db.verify_database(logger)
        metadata.clear()
        metadata.reflect()
        assert 'eph' in metadata.tables.keys()

        class Test(schema.Base):
            __tablename__ = 'test'
            id = sa.Column(sa.Integer, primary_key=True)

        assert 'test' not in metadata.tables.keys()
        db.verify_database(logger, names=['test'])
        metadata.clear()
        metadata.reflect()
        assert 'test' in metadata.tables.keys()

    def test_create_test_db(self, db):
        # just exercise the code
        assert db is not None

    def test_add_ephemeris_fixed(self, db):
        # already added via test database
        c = db.session.query(Eph).count()
        assert c == 3  # default number in test database

    def test_add_ephemeris_adaptable(self, db):
        c0 = db.session.query(Eph).count()
        db.add_ephemeris(2, '500', 2457799.5, 2457809.5, step=None,
                         source='mpc', cache=True)
        c1 = db.session.query(Eph).count()
        assert c1 - c0 == 36

    def test_add_found(self, db):
        rows = db.get_found().all()
        assert len(rows) == 3

        foundids = db.add_found_by_id(2, [3], '500', cache=True,
                                      update=False)
        assert len(foundids) == 0

        foundids = db.add_found_by_id(2, [3], '500', cache=True,
                                      update=True)
        assert len(foundids) == 1

    def test_add_found_by_id(self, db):
        # check test db found items
        rows = db.get_found().all()
        assert len(rows) == 3

        foundids = db.add_found_by_id(2, [3], '500', cache=True,
                                      update=False)
        assert len(foundids) == 0

        foundids = db.add_found_by_id(2, [3], '500', cache=True, update=True)
        assert len(foundids) == 1

    def test_add_object(self, db):
        # check test db object
        obj = db.get_objects().filter_by(desg="C/1995 O1").one()
        assert obj.objid == 1
        assert obj.desg == 'C/1995 O1'

    def test_add_object_error(self, db):
        with pytest.raises(ValueError):
            db.add_object(1)

    def test_add_observations(self, db):
        logger = Logger('test')

        c = db.session.query(Obs).count()
        assert c == N_tiles**2

        # add rows
        obs = [Obs(
            jd_start=2450005.0,
            jd_stop=2450010.0,
            fov='SRID=40001;POLYGON((0 0, 1 1, 1 -1, -1 -1, -1 1, 0 0))',
        )]
        db.add_observations(obs, logger=logger)
        c = db.session.query(Obs).count()
        assert c == N_tiles**2 + 1

    def test_add_observations_with_source(self, db):
        class Survey(Obs):
            __tablename__ = 'survey'
            id = sa.Column(sa.Integer, primary_key=True)
            obsid = sa.Column(sa.Integer, sa.ForeignKey('obs.obsid'))

        with db.session.connection() as con:
            schema.Base.metadata.create_all(con)

        fov = FieldOfView(RADec(np.random.rand(4, 2), unit='deg'))
        obs = [Survey(
            jd_start=2458300.5,
            jd_stop=2458300.51,
            fov=str(fov)
        )]
        db.add_observations(obs)

        c = db.session.query(Obs).count()
        assert c == N_tiles**2 + 1

        c = db.session.query(Survey).count()
        assert c == 1

    def test_clean_ephemeris(self, db):
        jda, jdb = 2458119.5, 2458121.5
        c = db.get_ephemeris(2, jda, jdb).count()
        assert c == 3
        c = db.clean_ephemeris(2, jda, jdb)
        assert c == 3
        c = db.get_ephemeris(2, jda, jdb).count()
        assert c == 0

    def test_clean_found(self, db):
        count = db.clean_found(2, None, None)
        assert count == 3

    # tested up to here
    def test_get_ephemeris_date_range(self, db):
        jd_range = db.get_ephemeris_date_range()
        assert np.allclose(jd_range, [2458119.5, 2458121.5])

        jd_range = db.get_ephemeris_date_range(objids=[2])
        assert np.allclose(jd_range, [2458119.5, 2458121.5])

    def test_get_ephemeris_date_range_error(self, db):
        with pytest.raises(NoEphemerisError):
            jd_range = db.get_ephemeris_date_range([3])

    def test_get_ephemeris(self, db):
        eph = db.get_ephemeris(2, 2458119.5, 2458121.5).all()
        assert len(eph) == 3

    def test_get_ephemeris_all(self, db):
        eph = db.get_ephemeris(None, None, None).all()
        assert len(list(eph)) == 3

    def test_get_ephemeris_interp(self, db):
        jdc = 2458120.0
        jda, jdb = 2458119.5, 2458120.5
        eph = db.get_ephemeris(2, jda, jdb).all()
        a = RADec(eph[0].ra, eph[0].dec, unit='rad')
        b = RADec(eph[1].ra, eph[1].dec, unit='rad')
        test = util.spherical_interpolation(a, b, jda, jdb, jdc)

        eph, v = db.get_ephemeris_interp(2, [jdc])
        assert np.isclose(eph.separation(test).value, 0)

    def test_get_found_by_id(self, db):
        c = db.get_found_by_id([1, 3]).count()
        assert c == 2

        assert db.get_found_by_id([100]).count() == 0

    def test_get_found_date(self, db):
        start = 2458119.5
        stop = start + 60 / 86400

        c = db.get_found(start=start, stop=stop).count()
        assert c == 1

    def test_get_found_by_obsid(self, db):
        c = db.get_found_by_obsid([1]).count()
        assert c == 0

    def test_get_found_by_object(self, db):
        c = db.get_found(obj=1).count()
        assert c == 0

        c = db.get_found(obj=2).count()
        assert c == 3

    def test_get_objects(self, db):
        rows = db.get_objects().all()
        objid = [r.objid for r in rows]
        desg = [r.desg for r in rows]
        assert len(objid) == 2
        assert 'C/1995 O1' in desg
        assert '2P' in desg

    def test_get_observation_date_range(self, db):
        jd_range = db.get_observation_date_range()
        assert np.allclose(jd_range, [2458119.5, 2458119.5 + 50 / 1440])

        fov = str(FieldOfView(RADec(np.random.rand(4, 2), unit='deg')))
        obs = []
        obs.append(Obs(jd_start=2458300.5, jd_stop=2458300.51, fov=fov))
        obs.append(Obs(jd_start=2458100.5, jd_stop=2458100.51, fov=fov))

        jd_range = db.get_observation_date_range(source='unspecified')
        assert np.allclose(jd_range, [2458119.5, 2458119.5 + 50 / 1440])

    def test_get_observation_date_range_error(self, db):
        with pytest.raises(SourceNotFoundError):
            jd_range = db.get_observation_date_range(source='blah')

        obs = db.get_observations_by_date().delete()
        with pytest.raises(SourceNotFoundError):
            jd_range = db.get_observation_date_range()

    def test_get_observations_by_id(self, db):
        c = db.get_observations_by_id([]).count()
        assert c == 0

        c = db.get_observations_by_id([1, 2, 3]).count()
        assert c == 3

        c = db.get_observations_by_id([100]).count()
        assert c == 0

    def test_get_observations_by_date(self, db):
        c = db.get_observations_by_date(2458119.5, 2458121.5).count()
        assert c == N_tiles**2

    def test_get_observations_covering(self, db):
        eph = db.get_ephemeris(2, None, None).all()
        coords = RADec.from_eph(eph)
        fov = str(FieldOfView(RADec.from_eph(eph)))

        epochs = [e.jd for e in eph]
        start = min(epochs)
        stop = max(epochs)

        nobs = (db.get_observations_covering(fov, start=start, stop=stop)
                .count())

        # for N_tiles == 10, ephemeris will fully be in just one box
        assert nobs == 1

    def test_get_observations_intersecting(self, db):
        eph = db.get_ephemeris(2, None, None).all()
        coords = RADec.from_eph(eph)
        fov = str(FieldOfView(RADec.from_eph(eph)))

        epochs = [e.jd for e in eph]
        start = min(epochs)
        stop = max(epochs)

        nobs = (db.get_observations_intersecting(fov, start=start, stop=stop)
                .count())

        # for N_tiles == 10, ephemeris will fully be in just one box
        assert nobs == 1

    def test_resolve_objects(self, db):
        objid, desg = list(zip(*db.resolve_objects([1, '2P'])))
        assert objid[0] == 1
        assert objid[1] == 2
        assert desg[0] == 'C/1995 O1'
        assert desg[1] == '2P'

    def test_resolve_object(self, db):
        objid, desg = db.resolve_object(1)
        assert objid == 1
        assert desg == 'C/1995 O1'

        objid, desg = db.resolve_object('2P')
        assert objid == 2
        assert desg == '2P'

        objid, desg = db.resolve_object('1P')
        assert objid is None

    def test_resolve_object_fail(self, db):
        with pytest.raises(BadObjectID):
            db.resolve_object(23)

max_stars_repo_name: Small-Bodies-Node/sbsearch

max_issues_repo_name: Small-Bodies-Node/sbsearch

max_forks_repo_name: Small-Bodies-Node/sbsearch

reason: import numpy,import astropy


text:
import numpy as np
from coptim.function import Function

T = np.array([15, 25, 35, 45, 55, 65, 75, 85, 105, 185, 245, 305, 365])
Y = np.array([0.038, 0.085, 0.1, 0.103, 0.093, 0.095, 0.088, 0.08, 0.073, 0.05, 0.038, 0.028, 0.02])

class Bateman(Function):

    def eval(self, x, t_list=T):
        y_x_t = x[2] * (np.exp(-x[0] * t_list) - np.exp(-x[1] * t_list))
        f_raw = np.power(y_x_t - Y, 2)

        return 0.5 * np.sum(f_raw)

    def gradient(self, x):
        y_x_t = x[2] * (np.exp(-x[0] * T) - np.exp(-x[1] * T))

        dx1 = 0
        dx2 = 0
        dx3 = 0
        for i in range(13):
            dx1 = dx1 + (x[2] * -T[i] * (np.exp(-x[0] * T[i]))) * (y_x_t[i] - Y[i])
            dx2 = dx2 + (x[2] * T[i] * np.exp(-x[1] * T[i])) * (y_x_t[i] - Y[i])
            dx3 = dx3 + (np.exp(-x[0] * T[i]) - np.exp(-x[1] * T[i])) * (y_x_t[i] - Y[i])

        return np.array([dx1, dx2, dx3])

    def hessian(self, x):
        pass

max_issues_repo_name: cmazzaanthony/coptim

max_forks_repo_name: cmazzaanthony/coptim