text:
#!/usr/bin/env python3

"""
Benchmark recovery of fuzzball contacts in the context of a match

"""

import os
import re
import math
import pickle
from pprint import pprint

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

import pyrosetta
from pyrosetta import rosetta

from ..utils import *
from ..design import generate_fuzzball_contact_rotamersets

def fuzzball_contact_recovery_bm(args):
    """
    Count recovered fuzzball contacts within the context of a match at different rotamer sampling levels

    Usage:
      bsff benchmark fuzzball_contact_recovery bulk <user_defined_dir> <match_dir> <fuzzball_dir> [options]
      bsff benchmark fuzzball_contact_recovery standalone <design_dir> <ligand_conformer_path> [options]

    Arguments:
      <user_defined_dir>        Path to project root directory
      <match_dir>               Directory containing matches to benchmark
      <fuzzball_dir>            Directory containing fuzzballs for current match set
      <ligand_conformer_path>   Path to ligand generated by molfile_to_params.py for current complex

    Options:
      -d, --dump_rotsets                Dump a PDB for the RotamerSets generated for each extrachi/sample level
      -s, --design_setup                Perform recovery on a complex where files are already setup for design
      -t=<taskid>, --taskid=<taskid>    Run benchmark for the <taskid> match in os.listdir(<match_dir>)
    """

    def perform_rotamer_recovery(match_pdb_path, design_setup=False, dump_rotamerset_pdb=False):

        # Set up infrastructure
        rotamer_info_list = list()
        rotamer_stats = dict()  # Dict containing stats for all rotamers evaluated for viable rotamers

        # Get matched motif residue and ligand positions
        match_dir, match_name = os.path.split(os.path.normpath(match_pdb_path))

        # Configure Paths for design setups
        if design_setup:
            params_dir = os.path.join(match_dir, 'Conformers')
            params_file = [file for file in os.listdir(params_dir) if file.endswith('.params')][0]  # Jank
            params_path = os.path.join(params_dir, params_file)
            fuzzball_index = 0
            iter_index = 0
            ligand_conformer_path = conformer_path
        else:
            conformer_name, fuzzball_resnums = find_conformer_and_constraint_resnums(match_name)
            params_path = os.path.join(user_defined_dir, 'Inputs', 'Rosetta_Inputs', f'{conformer_name}.params')
            ligand_conformer_path = os.path.join(user_defined_dir, 'Inputs', 'Rosetta_Inputs', f'{conformer_name}.pdb')

            # Get fuzzball index from match_pdb
            match_pdb_split = re.split('[-_]', match_name)
            fuzzball_index = match_pdb_split[match_pdb_split.index('fuzz') + 1]
            iter_index = match_pdb_split[match_pdb_split.index('iter') + 1]

        # Instantiate PyRosetta things
        my_options = [
            f"-extra_res_fa {params_path}",
            "-preserve_header",
            "-extrachi_cutoff=0",
            "-mute core.chemical.ResidueType core.pack.task"
        ]
        pyrosetta.init(options=' '.join(my_options))
        sfxn = rosetta.core.scoring.get_score_function()
        sfxn_weights = sfxn.weights()

        # Load match into pose
        match_pose = rosetta.core.pose.Pose()
        rosetta.core.import_pose.pose_from_file(match_pose, match_pdb_path)
        ligand_resnum = match_pose.size()

        # --- Load viable scaffold positions and corresponding residue types --- #

        match_residue_map_temp = pickle.load(open(os.path.join(fuzzball_dir, f'iter_{iter_index}-fuzz_{fuzzball_index}-match_residue_map.pickle'), 'rb'))

        # Only consider non-CGP positions
        match_residue_map = dict()
        for position, residues in match_residue_map_temp.items():
            if match_pose.residue(int(position)).name1() not in ['C', 'G', 'P']:
                match_residue_map[int(position)] = residues

        # --- Generate rotamers at different sampling levels --- #

        for extrachi in range(1, 5):

            # Set extrachi options
            for chi_on in range(1, extrachi + 1):
                rosetta.basic.options.set_boolean_option(f'packing:ex{chi_on}', True)

            # Track
            rotamer_stats[extrachi] = dict()

            for sample_level in range(1, 8):

                # Set sample level options for active chi
                for chi_sample in range(1, extrachi + 1):
                    rosetta.basic.options.set_integer_option(f'packing:ex{chi_sample}:level', sample_level)

                # Sanity check
                print('\n\nCURRENT OPTIONS:')
                for chi_check in range(1, 5):
                    print(f'packing:ex{chi_check}', rosetta.basic.options.get_boolean_option(f'packing:ex{chi_check}'))
                    print(f'packing:ex{chi_check}:level', rosetta.basic.options.get_integer_option(f'packing:ex{chi_check}:level'))
                print('\n\n')

                # Generate Rotamers
                viable_rotamers, rotamer_trial_info = generate_fuzzball_contact_rotamersets(ligand_conformer_path, match_pdb_path, match_pose, sfxn, match_residue_map,
                                                                                flag_special_rot=False, report_stats=True, RMSD_limit=1, dump_rotamerset_pdb=dump_rotamerset_pdb)
                rotamer_stats[extrachi][sample_level] = rotamer_trial_info

                # Match pose clone for rotamer scores
                match_pose_clone = match_pose.clone()
                sfxn(match_pose_clone)

                # Add information to a list of dicts (for making a dataframe later)
                for position in viable_rotamers:
                    for residue_type in viable_rotamers[position]:
                        for rotamer in viable_rotamers[position][residue_type]:

                            # Add rotamer to pose, score, get one body energies from energy graph
                            match_pose_clone.replace_residue(position, rotamer, False)
                            sfxn(match_pose_clone)
                            e_edge = match_pose_clone.energies().energy_graph()
                            current_edge = e_edge.find_energy_edge(position, ligand_resnum)

                            if current_edge is not None:
                                emap = current_edge.fill_energy_map()
                                bsff_score = sum([emap[rosetta.core.scoring.fa_rep] * 0.55,
                                                  emap[rosetta.core.scoring.fa_atr],
                                                  emap[rosetta.core.scoring.hbond_sc],
                                                  emap[rosetta.core.scoring.fa_sol],
                                                  emap[rosetta.core.scoring.fa_elec]
                                                  ])
                                total_score = current_edge.dot(sfxn_weights)
                            else:
                                print(f'No edge between {position} and {ligand_resnum}')
                                bsff_score = None
                                total_score = None

                            rotamer_fa_dun = float(match_pose_clone.energies().onebody_energies(position)[rosetta.core.scoring.fa_dun])

                            rotamer_data = {'extrachi': extrachi,
                                            'sample_level': sample_level,
                                            'position': position,
                                            'resname': residue_type,
                                            'fa_dun': rotamer_fa_dun,
                                            'bsff_score': bsff_score,
                                            'total_score': total_score,
                                            }

                            rotamer_info_list.append(rotamer_data)

        # if dump_rotamerset_pdb:
        #
        #     # --- Dump packer rotset (code adapted from design) --- #
        #
        #     # Ligand RS
        #     matched_ligand_rs = rosetta.core.select.residue_selector.ResidueIndexSelector(match_pose.size())
        #
        #     # Positions directly adjacent to ligand RS
        #     ligand_adjacent_rs = rosetta.core.select.residue_selector.NeighborhoodResidueSelector(matched_ligand_rs, 8, include_focus_in_subset=True)
        #
        #     # NATRO positions in pose RS
        #     natro_rs = rosetta.core.select.residue_selector.NotResidueSelector(ligand_adjacent_rs)
        #
        #     # Create and Populate Task Factory
        #
        #     task_factory = rosetta.core.pack.task.TaskFactory()
        #
        #     racaa = rosetta.core.pack.task.operation.RestrictAbsentCanonicalAASRLT()
        #     racaa.aas_to_keep('ADEFHIKLMNQRSTVWY')  # No CGP
        #     design_op = rosetta.core.pack.task.operation.OperateOnResidueSubset(racaa, ligand_adjacent_rs)
        #     task_factory.push_back(design_op)
        #
        #     natro_op = rosetta.core.pack.task.operation.OperateOnResidueSubset(
        #         rosetta.core.pack.task.operation.PreventRepackingRLT(), natro_rs)
        #     task_factory.push_back(natro_op)
        #
        #     fixed_ligand_op = rosetta.core.pack.task.operation.OperateOnResidueSubset(
        #         rosetta.core.pack.task.operation.PreventRepackingRLT(), matched_ligand_rs)
        #     task_factory.push_back(fixed_ligand_op)
        #
        #     # Extra rotamers
        #     extra_rotamers_op = rosetta.core.pack.task.operation.ExtraRotamersGeneric()
        #     extra_rotamers_op.ex1(True)
        #     extra_rotamers_op.ex2(True)
        #     extra_rotamers_op.ex1_sample_level(1)
        #     extra_rotamers_op.ex2_sample_level(1)
        #     task_factory.push_back(extra_rotamers_op)
        #
        #     # Create RotamerSets including fuzzball rotamers
        #
        #     design_packer_task = task_factory.create_task_and_apply_taskoperations(task_factory, match_pose)
        #
        #     sfxn.setup_for_packing(match_pose, design_packer_task.repacking_residues(), design_packer_task.designing_residues())
        #     packer_neighbor_graph = rosetta.core.pack.create_packer_graph(match_pose, sfxn, design_packer_task)
        #
        #     rotamer_sets = rosetta.core.pack.rotamer_set.RotamerSetsFactory.create_rotamer_sets(match_pose)
        #     rotamer_sets.set_task(design_packer_task)
        #     rotamer_sets.initialize_pose_for_rotsets_creation(match_pose)
        #     rotamer_sets.build_rotamers(match_pose, sfxn, packer_neighbor_graph)
        #
        #     rotamer_sets.dump_pdb(match_pose, f"{match_name.split('.')[0]}-packer_rotset.pdb")

        # Generate dataframe for vaiable rotamers at different sample levels
        df = pd.DataFrame(rotamer_info_list)
        df.to_csv(f"{match_name.split('.')[0]}-contact_recovery_benchmark.csv")

        # Keep track of sampling types, rotamers generated, and rotamers accepted
        with open(f"{match_name.split('.')[0]}-rotamer_statistics.pickle", 'wb') as rotamer_stats_pickle:
            pickle.dump(rotamer_stats, rotamer_stats_pickle)

    # Parse args
    user_defined_dir = args['<user_defined_dir>']
    match_directory = args['<design_dir>'] if args['<design_dir>'] else args['<match_dir>']
    fuzzball_dir = os.path.join(args['<design_dir>'], 'Fuzzballs') if args['<design_dir>'] else args['<fuzzball_dir>']
    dump_rotsets = args['--dump_rotsets']
    conformer_path = args['<ligand_conformer_path>']

    # Either run in series or in parallel
    if args['--taskid']:
        current_match = os.listdir(match_directory)[int(args['--taskid']) - 1]
        perform_rotamer_recovery(os.path.join(match_directory, current_match), design_setup=args['standalone'], dump_rotamerset_pdb=dump_rotsets)

    else:
        for match_pdb in pdb_check(match_directory):
            perform_rotamer_recovery(match_pdb, design_setup=args['standalone'], dump_rotamerset_pdb=dump_rotsets)

def rotamer_recovery_figures_for_tanja(args):
    """
    Yeah...

    Usage:
      bsff benchmark contact_recovery_figs <results_dir> [options]

    Arguments:
      <results_dir>           Directory containing csvs generated by fuzzball_contact_recovery_bm()

    Options:
      --skip_figures, -s      Skip making figures
      --help, -h              Show this message.

    """

    results_dir = args['<results_dir>']
    template_dir = os.path.join(os.path.dirname(__file__), '../..', 'Additional_Files', 'templates')

    from jinja2 import Environment, FileSystemLoader

    # http://eosrei.net/articles/2015/11/latex-templates-python-and-jinja2-generate-pdfs
    latex_jinja_env = Environment(
        block_start_string='\BLOCK{',
        block_end_string='}',
        variable_start_string='\VAR{',
        variable_end_string='}',
        comment_start_string='\#{',
        comment_end_string='}',
        line_statement_prefix='%%',
        line_comment_prefix='%#',
        trim_blocks=True,
        autoescape=False,
        loader=FileSystemLoader(template_dir))

    template = latex_jinja_env.get_template('rotamer_recovery_bm.tex')

    # Get unique matches in results dir
    match_set = set([match[:match.rfind('-')] for match in os.listdir(results_dir)])
    print(match_set)

    for unique_match in match_set:
        if not args['--skip_figures']:
            # todo: add option for output dir
            figure_output_dir = generate_figures(results_dir, unique_match)

    # Pass jinja dict with match_names as key and list of figures as value
    jijna_figure_dict = {match: dict() for match in match_set}
    for key, value in jijna_figure_dict.items():
        key_figure_list = [fig for fig in os.listdir(figure_output_dir) if key in fig]
        jijna_figure_dict[key] = {re.split('[-.]', fig)[-2]: os.path.join(results_dir, fig) for fig in key_figure_list}

    pprint(jijna_figure_dict)

    with open(f'Complementary_RotamerSet-analysis.tex', 'w') as asdf:
        asdf.write(template.render(jijna_figure_dict=jijna_figure_dict))

def generate_figures(results_dir, unique_match, figure_output_dir='benchmark_figures'):
    """
    Make pretty figures using data generated during a fuzzball assembly iteration.

    :param results_dir: Directory containing contact_recovery_benchmark.csv and rotamer_statistics.pickle for <unique_match>.
    :param unique_match: Name of the nucleated match scaffold.
    """

    os.makedirs(figure_output_dir, exist_ok=True)
    viable_residue_df = pd.read_csv(os.path.join(results_dir, f'{unique_match}-contact_recovery_benchmark.csv'))
    all_rotamers_stats = pickle.load(open(os.path.join(results_dir, f'{unique_match}-rotamer_statistics.pickle'), 'rb'))

    # --- BSFF vs Actual REU --- #
    """
    Seaborn lmplot for each extrachi and sample level
    """
    if viable_residue_df.empty is False:

        sns_reu_comparison = sns.FacetGrid(viable_residue_df, col="extrachi", row="sample_level")
        sns_reu_comparison = sns_reu_comparison.map(plt.scatter, "bsff_score", "total_score", edgecolor="w")
        sns_reu_comparison.savefig(os.path.join(figure_output_dir, f'{unique_match}-reu_comparison.png'), format='png', dpi=300, bbox_inches='tight')

    # --- rotamer recovery for sample level + extrachi --- #
    """
    Heatmap for the time being...
    """
    if viable_residue_df.empty is False:

        for extrachi in all_rotamers_stats:
            for sample_level in all_rotamers_stats[extrachi]:

                total_trial = 0
                total_accepted = 0

                for position in all_rotamers_stats[extrachi][sample_level]:
                    for contact_residue in all_rotamers_stats[extrachi][sample_level][position]:
                        total_trial += all_rotamers_stats[extrachi][sample_level][position][contact_residue]['num_rotamers']
                        total_accepted += all_rotamers_stats[extrachi][sample_level][position][contact_residue][
                            'rotamers_accepted']

        rot_recovery_fig, rot_recovery_ax = plt.subplots()

        # Aggregate counts
        heatmap_counts = viable_residue_df.groupby(['extrachi', 'sample_level']).agg('count').pivot_table(index='extrachi',
                                                                                columns='sample_level',
                                                                                values='total_score')
        print(heatmap_counts)

        sns.heatmap(heatmap_counts,
                    cmap="Blues", cbar=False,
                    # norm=matplotlib.colors.LogNorm(vmin=heatmap_counts.min().min(), vmax=heatmap_counts.max().max()),
                    vmin=heatmap_counts.min().min(),
                    vmax=heatmap_counts.max().max(),
                    square=True,
                    annot=True, fmt='g', annot_kws={'color': 'black', 'fontsize': 'small'},
                    ax=rot_recovery_ax)
        rot_recovery_fig.savefig(os.path.join(figure_output_dir, f'{unique_match}-rot_recovery.png'), format='png', dpi=300, bbox_inches='tight')

    # --- fa_dun vs sample level + extrachi --- #
    """
    Look at fa_dun as extrachi/sampling increases
    """
    if viable_residue_df.empty is False:

        fa_dun_fig, fa_dun_ax = plt.subplots(figsize=(4, 12))
        viable_residue_df['extrachi_sample'] = [f'Extrachi {extrachi}, Sampling {sample}' for extrachi, sample in
                                                zip(viable_residue_df['extrachi'], viable_residue_df['sample_level'])]
        sns.violinplot(x='fa_dun', y='extrachi_sample', data=viable_residue_df, scale='width', inner='quartile',
                       ax=fa_dun_ax)
        fa_dun_fig.savefig(os.path.join(figure_output_dir, f'{unique_match}-fa_dun.png'), format='png', dpi=300, bbox_inches='tight')

    # --- Energy distribution of accepted residues --- #

    if viable_residue_df.empty is False:

        reu_fig, reu_ax = plt.subplots(figsize=(4, 12))
        viable_residue_df['extrachi_sample'] = [f'Extrachi {extrachi}, Sampling {sample}' for extrachi, sample in
                                                zip(viable_residue_df['extrachi'], viable_residue_df['sample_level'])]
        sns.violinplot(x='total_score', y='extrachi_sample', data=viable_residue_df, scale='width', inner='quartile',
                       ax=reu_ax)
        reu_fig.savefig(os.path.join(figure_output_dir, f'{unique_match}-reu_accepted_dist.png'), format='png', dpi=300, bbox_inches='tight')

    # --- RMSD vs energy (highlight Accepted/Rejected) --- #

    # Construct dataframe fro pickle dict
    df_dict_list = list()

    for extrachi in all_rotamers_stats:
        for sample_level in all_rotamers_stats[extrachi]:
            for position in all_rotamers_stats[extrachi][sample_level]:
                for contact_residue in all_rotamers_stats[extrachi][sample_level][position]:
                    for rotamer_info in all_rotamers_stats[extrachi][sample_level][position][contact_residue]['rotamer_info']:
                        row_dict = {'extrachi': extrachi,
                                    'sample_level': sample_level,
                                    'position': position,
                                    'contact_residue': contact_residue,
                                    'min_RMSD': min(rotamer_info[0]),
                                    'contact_REU': rotamer_info[2],
                                    'accepted': True if (min(rotamer_info[0]) < 1 and rotamer_info[2] < 5) else False,
                                    }

                        df_dict_list.append(row_dict)

    residue_rmsd_reu_df = pd.DataFrame(df_dict_list)
    select_accepted_group = True in residue_rmsd_reu_df['accepted']

    if viable_residue_df.empty is False:
        rmsd_reu_grid = sns.FacetGrid(residue_rmsd_reu_df.groupby('accepted').get_group(select_accepted_group), col="extrachi", row="sample_level")
        rmsd_reu_grid = rmsd_reu_grid.map(plt.scatter, "min_RMSD", "contact_REU", edgecolor="w")
    else:
        rmsd_reu_grid = sns.FacetGrid(residue_rmsd_reu_df, col="extrachi", row="sample_level", sharey=False)
        rmsd_reu_grid = rmsd_reu_grid.map(plt.hist, "min_RMSD", edgecolor="w")

    rmsd_reu_grid.savefig(os.path.join(figure_output_dir, f'{unique_match}-rmsd_reu_grid.png'), format='png', dpi=300, bbox_inches='tight')

    # --- FacetGrid Accepted/Rejected at each position --- #
    accepted_grid = sns.FacetGrid(residue_rmsd_reu_df, col="extrachi", row="sample_level", hue='accepted')
    accepted_grid = accepted_grid.map(plt.hist, "position", bins=list(residue_rmsd_reu_df['position'].unique()))
    accepted_grid.set(yticks=[10 ** i for i in range(1, 6)], yticklabels=[10 ** i for i in range(1, 6)])

    for x_idx, x_ax in enumerate(accepted_grid.axes):
        for y_idx, y_ax in enumerate(accepted_grid.axes[x_idx]):
            accepted_grid.axes[x_idx][y_idx].set_yscale('log')

    accepted_grid.add_legend()
    accepted_grid.savefig(os.path.join(figure_output_dir, f'{unique_match}-accepted_grid.png'), format='png', dpi=300, bbox_inches='tight')

    # --- How many contact modes are being fulfilled? --- #
    """
    I'll do this later... maybe...
    """

    return figure_output_dir

max_stars_repo_name: jaaamessszzz/BindingSitesFromFragments

max_issues_repo_name: jaaamessszzz/BindingSitesFromFragments

max_forks_repo_name: jaaamessszzz/BindingSitesFromFragments


text:
# 求解定积分 x^2 区间[1, 2]; 投点法
import numpy as np
import tqdm
sum = 0
for i in tqdm.tqdm(range(3000000)):

        random_generate = np.array([np.random.uniform(1, 2), np.random.uniform(0, 4)])
        if np.square(random_generate[0]) > random_generate[1]:
                sum += 1
print(sum)
area = 4 * sum / 3000000
print('Area is:{}'.format(area))

ext: py

reason: import numpy


text:
#!/usr/bin/env python3
import unittest

import tensorflow as tf
import termcolor
import numpy as np
from numpy import exp, log
from numpy.testing import (assert_allclose, assert_approx_equal,
                           assert_array_equal, assert_raises)

from reward_predictor import RewardPredictorNetwork
from reward_predictor_core_network import net_cnn

class TestRewardPredictor(unittest.TestCase):

    def setUp(self):
        self.create_reward_predictor(dropout=0.5, batchnorm=True)
        termcolor.cprint(self._testMethodName, 'red')

    def create_reward_predictor(self, dropout, batchnorm):
        tf.reset_default_graph()
        self.sess = tf.Session()
        self.rpn = RewardPredictorNetwork(batchnorm=batchnorm, dropout=dropout,
                                          lr=1e-3,
                                          core_network=net_cnn)
        self.sess.run(tf.global_variables_initializer())

    def test_weight_sharing(self):
        """
        Check that both legs of the network give the same reward output
        for the same segment input.
        """

        # We deliberately /don't/ use the same dropout for each leg of the
        # network. (If we do use the same dropout, without batchnorm,
        # Pong doesn't train successfully. If we use different dropout, Pong
        # does train successfully. I haven't tried training Pong with
        # batchnorm.) So we disable dropout for this test.
        self.create_reward_predictor(dropout=0.0, batchnorm=True)

        s = 255 * np.random.rand(100, 84, 84, 4)
        feed_dict_nontraining = {
            self.rpn.s1: [s],
            self.rpn.s2: [s],
            self.rpn.training: True
        }
        feed_dict_training = {
            self.rpn.s1: [s],
            self.rpn.s2: [s],
            self.rpn.training: False
        }
        for feed_dict in [feed_dict_nontraining, feed_dict_training]:
            for _ in range(3):  # to check different dropouts
                [rs1], [rs2] = self.sess.run([self.rpn.rs1, self.rpn.rs2], feed_dict)
                # Check rs1 != 0.0
                assert_raises(AssertionError, assert_array_equal, rs1, 0.0)
                assert_allclose(rs1, rs2)

    def test_batchnorm_sharing(self):
        """
        Check that batchnorm statistics are the same between the two legs of
        the network.
        """
        n_frames = 20
        s1 = 255 * np.random.normal(loc=1.0, size=(n_frames, 84, 84, 4))
        s2 = 255 * np.random.normal(loc=-1.0, size=(n_frames, 84, 84, 4))
        feed_dict = {
            self.rpn.s1: [s1],
            self.rpn.s2: [s2],
            self.rpn.pref: [[0.0, 1.0]],
            self.rpn.training: True}
        self.sess.run(self.rpn.train, feed_dict)

        feed_dict = {self.rpn.s1: [s1], self.rpn.s2: [s1], self.rpn.training: False}
        [rs1], [rs2] = self.sess.run([self.rpn.rs1, self.rpn.rs2], feed_dict)
        # Check rs1 != 0.0
        assert_raises(AssertionError, assert_array_equal, rs1, 0.0)
        assert_allclose(rs1, rs2)

    def test_loss(self):
        """
        Check that the loss is calculated correctly.
        """
        # hack to ensure numerical stability
        rs1 = rs2 = 100
        n_frames = 20
        while rs1 > 50 or rs2 > 50:
            s1 = 255 * np.random.normal(loc=1.0, size=(n_frames, 84, 84, 4))
            s2 = 255 * np.random.normal(loc=-1.0, size=(n_frames, 84, 84, 4))
            feed_dict = {
                self.rpn.s1: [s1],
                self.rpn.s2: [s2],
                self.rpn.training: True
            }
            [rs1], [rs2] = self.sess.run([self.rpn.rs1, self.rpn.rs2],
                                         feed_dict)

        prefs = [[0.0, 1.0], [1.0, 0.0], [0.5, 0.5]]
        for pref in prefs:
            feed_dict[self.rpn.pref] = [pref]
            [rs1], [rs2], loss = self.sess.run(
                [self.rpn.rs1, self.rpn.rs2, self.rpn.loss], feed_dict)

            p_s1_s2 = exp(rs1) / (exp(rs1) + exp(rs2))
            p_s2_s1 = exp(rs2) / (exp(rs1) + exp(rs2))

            expected = -(pref[0] * log(p_s1_s2) + pref[1] * log(p_s2_s1))
            assert_approx_equal(loss, expected, significant=3)

    def test_batches(self):
        """
        Present a batch of two trajectories and check that we get the same
        results as if we'd presented the trajectories individually.
        """
        n_segs = 2
        n_frames = 20
        prefs = [[0., 1.], [1., 0.]]
        s1s = []
        s2s = []
        for _ in range(n_segs):
            s1 = 255 * np.random.normal(loc=1.0, size=(n_frames, 84, 84, 4))
            s2 = 255 * np.random.normal(loc=-1.0, size=(n_frames, 84, 84, 4))
            s1s.append(s1)
            s2s.append(s2)

        # Step 1: present all trajectories as one big batch
        feed_dict = {
            self.rpn.s1: s1s,
            self.rpn.s2: s2s,
            self.rpn.pref: prefs,
            self.rpn.training: False
        }
        rs1_batch, rs2_batch, pred_batch, loss_batch = self.sess.run(
            [self.rpn.rs1, self.rpn.rs2, self.rpn.pred, self.rpn.loss],
            feed_dict)

        # Step 2: present trajectories individually
        rs1_nobatch = []
        rs2_nobatch = []
        pred_nobatch = []
        loss_nobatch = 0
        for i in range(n_segs):
            feed_dict = {
                self.rpn.s1: [s1s[i]],
                self.rpn.s2: [s2s[i]],
                self.rpn.pref: [prefs[i]],
                self.rpn.training: False
            }
            [rs1], [rs2], [pred], loss = self.sess.run(
                [self.rpn.rs1, self.rpn.rs2, self.rpn.pred, self.rpn.loss],
                feed_dict)
            rs1_nobatch.append(rs1)
            rs2_nobatch.append(rs2)
            pred_nobatch.append(pred)
            loss_nobatch += loss

        # Compare
        assert_allclose(rs1_batch, rs1_nobatch, atol=1e-5)
        assert_allclose(rs2_batch, rs2_nobatch, atol=1e-5)
        assert_allclose(pred_batch, pred_nobatch, atol=1e-5)
        assert_approx_equal(loss_batch, loss_nobatch, significant=4)

    def test_training(self):
        """
        Present two trajectories with different preferences and see whether
        training really does work (whether the reward predicted by the network
        matches the preferences after a few loops of running the training
        operation).

        Note: because of variations in training, this test does not always pass.
        """
        n_frames = 20
        s1 = 255 * np.random.normal(loc=1.0, size=(n_frames, 84, 84, 4))
        s2 = 255 * np.random.normal(loc=-1.0, size=(n_frames, 84, 84, 4))

        feed_dict = {
            self.rpn.s1: [s1],
            self.rpn.s2: [s2]
        }

        prefs = [[0.0, 1.0], [1.0, 0.0], [0.5, 0.5]]
        for pref in prefs:
            print("Preference", pref)
            feed_dict[self.rpn.pref] = [pref]
            # Important to reset batch normalization statistics
            self.sess.run(tf.global_variables_initializer())
            for _ in range(150):
                feed_dict[self.rpn.training] = True
                self.sess.run(self.rpn.train, feed_dict)
                # Uncomment these for more thorough manual testing.
                # (For the first case, rs1 should become higher
                #  than rs2, and the distance between them should increase;
                #  for the second case, rs2 should become higher;
                #  for the third case, they should become approximately the
                #  same.)
            """
                feed_dict[self.rpn.training] = False
                ops = [self.rpn.rs1, self.rpn.rs2, self.rpn.loss]
                [rs1], [rs2], loss = self.sess.run(ops, feed_dict)
                print(" ".join(3 * ["{:>8.3f}"]).format(rs1, rs2, loss))
            print()
            """

            feed_dict[self.rpn.training] = False
            [rs1], [rs2] = self.sess.run([self.rpn.rs1, self.rpn.rs2], feed_dict)

            if pref[0] > pref[1]:
                self.assertGreater(rs1 - rs2, 10)
            elif pref[1] > pref[0]:
                self.assertGreater(rs2 - rs1, 10)
            elif pref[0] == pref[1]:
                self.assertLess(abs(rs2 - rs1), 2)

    def test_training_batches(self):
        """
        Check that after training with a batch of 4 segments, each with their own preferences,
        the predicted preference for each of the segments is as expected.
        """
        n_frames = 20
        s1s = 255 * np.random.normal(loc=1.0, size=(4, n_frames, 84, 84, 4))
        s2s = 255 * np.random.normal(loc=-1.0, size=(4, n_frames, 84, 84, 4))
        prefs = [[1.0, 0.0], [1.0, 0.0], [0.0, 1.0], [0.0, 1.0]]
        feed_dict = {
            self.rpn.s1: s1s,
            self.rpn.s2: s2s,
            self.rpn.pref: prefs,
            self.rpn.training: True
        }

        for i in range(100):
            if i % 10 == 0:
                print("Training {}/100".format(i))
            self.sess.run(self.rpn.train, feed_dict)

        feed_dict[self.rpn.training] = False
        preds = self.sess.run(self.rpn.pred, feed_dict)
        assert_allclose(preds[0], [1., 0.], atol=1e-1)
        assert_allclose(preds[1], [1., 0.], atol=1e-1)
        assert_allclose(preds[2], [0., 1.], atol=1e-1)
        assert_allclose(preds[3], [0., 1.], atol=1e-1)

    def test_accuracy(self):
        """
        Test accuracy op.
        """
        n_frames = 20
        batch_n = 16
        s1s = 255 * np.random.normal(loc=1.0, size=(batch_n, n_frames, 84, 84, 4))
        s2s = 255 * np.random.normal(loc=-1.0, size=(batch_n, n_frames, 84, 84, 4))
        possible_prefs = [[1.0, 0.0], [0.0, 1.0]]
        possible_prefs = np.array(possible_prefs)
        prefs = possible_prefs[np.random.choice([0, 1], size=batch_n)]

        feed_dict = {
            self.rpn.s1: s1s,
            self.rpn.s2: s2s,
            self.rpn.pref: prefs,
            self.rpn.training: True
        }

        # Steer away from chance performance
        for _ in range(5):
            self.sess.run(self.rpn.train, feed_dict)

        feed_dict[self.rpn.training] = False
        preds = self.sess.run(self.rpn.pred, feed_dict)
        n_correct = 0
        for pref, pred in zip(prefs, preds):
            if pref[0] == 1.0 and pred[0] > pred[1] or \
               pref[1] == 1.0 and pred[1] > pred[0]:
                n_correct += 1
        accuracy_expected = n_correct / batch_n

        accuracy_actual = self.sess.run(self.rpn.accuracy, feed_dict)

        assert_approx_equal(accuracy_actual, accuracy_expected)

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

max_stars_repo_name: unghee/learning-from-human-preferences

max_issues_repo_name: unghee/learning-from-human-preferences

max_forks_repo_name: unghee/learning-from-human-preferences

reason: import numpy,from numpy


text:
from codemaker.embedding import SDAEmbedder
from codemaker.transforms import random_project
from codemaker.evaluation import Neighbors, local_match
import numpy as np
from nose.tools import *
from scikits.learn import datasets
import random

def test_compute_embedding(check_asserts=True):
    np.random.seed(0)
    random.seed(0)

    # sample data from the digits 8x8 pixels dataset
    digits = datasets.load_digits()
    data = digits.data
    n_samples, n_features = data.shape
    low_dim = 2

    # baseline score using a random 2D projection
    projection = random_project(data, target_dim=2, rng=np.random)
    score = local_match(data, projection, query_size=50, ratio=0.1, seed=0)
    assert_almost_equal(score, 0.12, 2)

    # compute an embedding of the data
    embedder = SDAEmbedder((n_features, 40, 15, low_dim),
                           noise=0.1,
                           reconstruction_penalty=.0,
                           embedding_penalty=1,
                           sparsity_penalty=0.0, learning_rate=0.1, seed=0)
    embedder.pre_train(data, epochs=500, batch_size=5)

    code = embedder.encode(data)
    assert_equal(code.shape, (n_samples, low_dim))

    # compare nearest neighbors
    score = local_match(data, code, query_size=50, ratio=0.1, seed=0)

    assert_almost_equal(score, 0.33, 1)

max_stars_repo_name: ogrisel/codemaker

max_issues_repo_name: ogrisel/codemaker

max_forks_repo_name: ogrisel/codemaker


text:
# Lint as: python3
"""Statistical models based on Neural Networks."""
from typing import Any, Callable
import dataclasses

from epi_forecast_stat_mech.statistical_models import base
from epi_forecast_stat_mech.statistical_models import probability as stat_prob
from epi_forecast_stat_mech.statistical_models import tree_util

from flax import nn
import jax
import jax.numpy as jnp

import tensorflow_probability as tfp
tfp = tfp.experimental.substrates.jax
tfd = tfp.distributions

class PerceptronModule(nn.Module):
  """Multi-layer perceptron network module."""

  def apply(
      self,
      inputs,
      output_size,
      hidden_layer_sizes=(64, 64),
      activation=jax.nn.relu,
      batch_norm_style='once'
  ):
    """Computes the output of a multi-layer perceptron give `inputs`."""
    x = inputs
    if batch_norm_style == 'none' or batch_norm_style == 'once':
      if batch_norm_style == 'once':
        x = nn.BatchNorm(x)
      for size in hidden_layer_sizes:
        x = nn.Dense(x, size)
        x = activation(x)
      return nn.Dense(x, output_size)
    elif batch_norm_style == 'layerwise_whiten':
      for size in hidden_layer_sizes:
        x = nn.BatchNorm(x, bias=False, scale=False)
        x = nn.Dense(x, size)
        x = activation(x)
      x = nn.BatchNorm(x, bias=False, scale=False)
      return nn.Dense(x, output_size)
    else:
      raise ValueError(f'Unexpected batch_norm_style: {batch_norm_style}')

class PlainLinearModule(
    PerceptronModule.partial(
        hidden_layer_sizes=(), activation=None, batch_norm_style='none')):
  pass

class LinearModule(
    PerceptronModule.partial(
        hidden_layer_sizes=(), activation=None, batch_norm_style='once')):
  pass

@dataclasses.dataclass
class NormalDistributionModel(base.StatisticalModel):
  """Statistical model based on normal distribution.

  Attributes:
    predict_module: flax.nn.Module that takes `inputs` and `output_size`
      arguments and returns array of shape `[batch, output_size]` that will be
      used to predict locations of the gaussian distributed `observations` and
      possibly scales, depending on whether `error_model` is 'full'.
    error_model: (string) Either 'full' to represent using the predict_module
      to estimate the error-scale (heteroscedastic). Or 'plugin' to
      represent that the scale should be based on a homoscedastic plugin.
    log_prior_fn: function that computes log_prior on a parameters of the
      `predict_module`.
    scale_eps: minimal value for the predicted scale if scale.
  """
  predict_module: nn.Module = PerceptronModule
  log_prior_fn: Callable[..., Any] = None
  error_model: str = 'full'
  scale_eps: float = 1e-2

  def _output_size_and_unpack_fn(self, output_structure):
    output_array, unpack_fn = tree_util.pack(output_structure)
    output_size = output_array.shape[-1]
    if self.error_model == 'full':
      output_size *= 2
    return output_size, unpack_fn

  def init_parameters(self, rng, inputs, output_structure):
    """Returns initial parameters generated at model construction time."""
    output_size, _ = self._output_size_and_unpack_fn(output_structure)
    return self.predict_module.init(rng, inputs, output_size)[1]

  def log_prior(self, parameters):
    """Returns the log probability of `parameters` based on priors.

    Args:
      parameters: parameters of the statistical model.

    Returns:
      log-probabilities for parameter.
    """
    if self.log_prior_fn is None:
      return 0.
    return self.log_prior_fn(parameters)

  def log_likelihood(self, parameters, covariates, observations):
    """Returns the log likelihood of `observations`.

    Args:
      parameters: parameters of the statistical model.
      covariates: A numpy array of shape "location" x "static_covariate".
      observations: A tree of mech_params that we want to explain.

    Returns:
      A log-likelihood. If self.error_model == 'full', the return is a tree of
      the shape of observations (the log_likelihood of each observation). In the
      'plugin' case, the return is a scalar (the sum).
    """
    if self.error_model == 'full':
      posterior = self.predict(parameters, covariates, observations)
      res = tree_util.tree_multimap(
          lambda p, o: p.log_prob(o), posterior, observations)
      return res
    elif self.error_model == 'plugin':
      loc, _ = self.get_loc_scale(parameters, covariates, observations)
      tree_error = tree_util.tree_multimap(
          lambda hat, o: o - hat, loc, observations)
      error, _ = tree_util.pack(tree_error)
      logprob = stat_prob.gaussian_error_logprob_with_bottom_scale(
          error, self.scale_eps, axis=0)
      return logprob
    else:
      raise ValueError(f'unexpected error_model: {self.error_model}')

  def get_loc_scale(self, parameters, covariates, observations):
    """Computes loc and scale for `observations` based on `covariates`.

    Args:
      parameters: parameters of the statistical model.
      covariates: array representing covariates for each location.
      observations: structure of observations to be predicted.

    Returns:
      pytree of probability distributions over `observations` given the
      `parameters` of the statistical model and `covariates`.
    """
    output_size, unpack_fn = self._output_size_and_unpack_fn(observations)
    raw_predictions = self.predict_module.call(
        parameters, covariates, output_size)
    if self.error_model == 'full':
      loc, raw_scale = jnp.split(raw_predictions, 2, -1)
      scale = jax.nn.softplus(raw_scale) + self.scale_eps
      return unpack_fn(loc), unpack_fn(scale)
    elif self.error_model == 'plugin':
      loc = raw_predictions
      return unpack_fn(loc), None
    else:
      raise ValueError(f'unexpected error_model: {self.error_model}')

  def predict(self, parameters, covariates, observations):
    """Predicts a distribution over `observations` based on `covariates`.

    Args:
      parameters: parameters of the statistical model.
      covariates: array representing covariates for each location.
      observations: structure of observations to be predicted.

    Returns:
      pytree of probability distributions over `observations` given the
      `parameters` of the statistical model and `covariates`.
    """
    loc, scale = self.get_loc_scale(parameters, covariates, observations)
    return tree_util.tree_multimap(
        lambda l, s: tfd.Normal(loc=l, scale=s),
        loc, scale)

  def linear_coefficients(self, parameters):
    dense_name = [x for x in parameters.keys() if 'Dense' in x][0]
    kernel = parameters[dense_name]['kernel']
    bias = parameters[dense_name]['bias']
    if self.error_model == 'full':
      kernel, _ = jnp.split(kernel, 2, -1)
      bias, _ = jnp.split(bias, 2, -1)
    return kernel, bias

max_stars_repo_name: HopkinsIDD/EpiForecastStatMech

max_issues_repo_name: HopkinsIDD/EpiForecastStatMech

max_forks_repo_name: HopkinsIDD/EpiForecastStatMech

reason: import jax


text:
from collections import Counter
import torch
import torch.nn as nn
import torch.nn.functional as F
import torch.optim as optim
import numpy as np

MAX_WORD_SIZE = 60000
BATCH_SIZE = 64
device = torch.device('cuda' if torch.cuda.is_available() else 'cpu')

def load_data(filename):
    x, y = [], []
    with open(filename, 'r', encoding='utf-8') as f:
        for line in f:
            line = line.strip().split('\t')
            x.append(line[1])
            y.append(line[0])
    return x, y

train_input, train_output = load_data('train.txt')
test_input, test_output = load_data('test.txt')

print(f'输入:{train_input[3][:20]},输出:{train_output[3]}')
print(f'输入:{test_input[3][:20]},输出:{test_output[3]}')

def build_vocab(text):
    all_c = ' '.join(text)
    all_c = all_c.split(' ')
    dic = Counter(all_c)
    dic = dic.most_common(MAX_WORD_SIZE-2)
    dic, _ = zip(*dic)
    index2word = ['<PAD>']+['<UNK>']+list(dic)
    word2index = {word: index for index, word in enumerate(index2word)}
    return index2word, word2index

index2word, word2index = build_vocab(train_input)
count_out = Counter(train_output)
index2label = list(count_out)
label2index = {w: i for i, w in enumerate(index2label)}

def encode(train_input, train_output, word2index, label2index):

    train_input1 = [[word2index.get(w, 1)
                     for w in i.split(' ')] for i in train_input]
    train_output1 = [[label2index.get(w, 1)
                      for w in i.split(' ')] for i in train_output]

    sorted_index = sorted(range(len(train_input1)),
                          key=lambda x: len(train_input1[x]))
    out_input = [train_input1[i] for i in sorted_index]
    out_output = [train_output1[i] for i in sorted_index]
    return out_input, out_output

train_in, train_out = encode(
    train_input, train_output, word2index, label2index)
print(' '.join([index2word[i] for i in train_in[1000]]))
print(' '.join([index2label[i] for i in train_out[1000]]))

test_in,test_out = encode(test_input,test_output,word2index,label2index)

def get_batch_index(n, batch_index_size):
    idx_list = np.arange(0, n, batch_index_size)
    np.random.shuffle(idx_list)
    batch_index = []
    for i in idx_list:
        batch_index.append(np.arange(i, min(batch_index_size+i, n)).tolist())
    return batch_index

def process_data(sentences, label):
    lens = [len(s) for s in sentences]
    n = len(sentences)
    max_len = np.max(lens)

    out = np.zeros((n, max_len)).astype('int32')
    for i in range(n):
        out[i, :lens[i]] = sentences[i]
    return out, np.array(lens).astype('int32'), np.array(label).astype('int32')

def get_batch(train_in, train_out, batch_size):
    batches_index = get_batch_index(len(train_in), batch_size)
    batches_data = []
    for batch in batches_index:
        sentences = [train_in[i]for i in batch]
        labal = [train_out[i][0] for i in batch]
        batches_data.append(process_data(sentences, labal))
    return batches_data

# [(train_in[i],train_out[i]) for i in range(3)]
train_data = get_batch(train_in, train_out, BATCH_SIZE)
# batch，句子，句子长度，标签
test_data = get_batch(test_in,test_out,BATCH_SIZE)
len(test_data[0])

class avgModule(nn.Module):
    def __init__(self, vocab_size, embedding_size, pad_idx, output_size):
        super(avgModule, self).__init__()
        self.embedding = nn.Embedding(vocab_size, embedding_size, pad_idx)
        self.fc = nn.Linear(embedding_size, output_size)
        self.dropout = nn.Dropout(0.5)

    def forward(self, input):
        # input : batch_size,vocab_size

        embeded = self.embedding(input)  # b_size,vocab_size,embeding_size
        embeded = self.dropout(embeded)
        p1 = F.avg_pool2d(embeded, (len(embeded[1]), 1)).squeeze()
        p2 = F.max_pool2d(embeded, (len(embeded[1]), 1)).squeeze()
        pooled=p1+p2
        return self.fc(pooled)

vocab_size = MAX_WORD_SIZE
embedding_size = 100
pad_idx = 0
output_size = 5
model = avgModule(vocab_size, embedding_size, pad_idx, output_size).to(device)

criterion = nn.CrossEntropyLoss().to(device)

optimizer = optim.Adam(model.parameters())

train_data[0][2]

def accuracy(preds, y):
    f_preds = preds.max(1)[1]
    correct = (f_preds == y).float()
    acc = sum(correct)/len(correct)
    return acc

from tensorboardX import SummaryWriter

writer = SummaryWriter('runs/exp')
def evaluate(model,criterion,data,epoch):
    epoch_acc = 0.0
    epoch_loss = 0.0
    xx =  -1
    for (x,length,y) in data:
        x = torch.from_numpy(x).long().to(device)
        y = torch.from_numpy(y).long().to(device)
        print(x.shape)
        preds = model(x)
        # if(xx==-1):
        #     writer.add_graph(model,input_to_model=x)
        #     xx=1
        loss = criterion(preds,y)
        acc = accuracy(preds,y)
        epoch_acc += acc
        epoch_loss += loss
    writer.add_scalar('loss',epoch_loss/len(data),epoch)
    writer.add_scalar('acc',epoch_acc/len(data),epoch)
    print(f'测试集：准确率:{epoch_acc/len(data)},loss:{epoch_loss/len(data)}')

def train(model, criterion, optimizer, data):
    epochs = 100
    for epoch in range(epochs):
        epoch_loss = 0
        epoch_acc = 0
        for (x, length, y) in data:
            x = torch.from_numpy(x).long().to(device)
            y = torch.from_numpy(y).long().to(device)
            preds = model(x)

            loss = criterion(preds, y)
            acc = accuracy(preds, y)

            optimizer.zero_grad()
            loss.backward()
            optimizer.step()
            epoch_loss += loss
            epoch_acc += acc
        with torch.no_grad():
            model.eval()
            evaluate(model,criterion,test_data,epoch)
            model.train()
        print(f'Epoch:{epoch},精准度:{epoch_acc/len(data)},loss:{epoch_loss/len(data)}')

train(model, criterion, optimizer, train_data)

evaluate(model,criterion,test_data,2)

[index2word[i] for i in test_data[0][0][1]]
tensor1 = torch.tensor(test_data[0][0][1]).long().cuda()
tensor1 = tensor1.view(1,-1)
tensor2 =  torch.tensor(test_data[0][0]).long().to(device)
pred = model(tensor2)
a = [index2label[i] for i in pred.max(1)[1]]
b = [index2label[i] for i in test_data[0][2]]
for i in zip(a,b):
    print(i)

t1 = torch.randn(2,3,4)
t1
F.max_pool2d(t1,(3,1))

t2 = torch.ones(3,4)
t1+t2

index2word[100:800]
it = iter(model.parameters())
i = next(it)

writer.add_embedding(mat=i[100:800],metadata=index2word[100:800],global_step=1)

max_stars_repo_name: aachong/SUDA-newStudents

max_issues_repo_name: aachong/SUDA-newStudents

max_forks_repo_name: aachong/SUDA-newStudents


text:
"""
A number of image preprocessing functions for the style transfer code
Author: Christian Baumgartner (c.f.baumgartner@gmail.com)
Date: 1. April 2016
"""

import numpy as np
import skimage.transform
from lasagne.utils import floatX

def convert_to_uint8(image):
    """
    Scales an image into the [0,255] range and casts it to uint8
    """
    image = image - image.min()
    image = 255.0*np.divide(image.astype(np.float32),image.max())
    return image.astype(np.uint8)

def check_image_format(image):
    """
    Checks and fixes for the following:
     - is image really uint8
     - does the image have an alpha channel
     - is the image gray scale
    """

    # pngs can sometimes not be uint8, so make sure
    image = convert_to_uint8(image)

    # delete alpha channel if there is one
    image = image[:,:,0:3]

    # add third axis to gray scale imaes
    if len(image.shape) == 2:
        image = image[:, :, np.newaxis]
        image = np.repeat(image, 3, axis=2)

    return image

def resize_image_to_input_size(image, im_size):
    """
    Resizes an images such that the smaller dimension equals im_size and the larger dimension
    is scaled keeping the aspect ration. im_size is given by the user as global parameter.
    """
    # resize image
    h, w, _ = image.shape
    if h < w:
        image = skimage.transform.resize(image, (im_size, w * im_size / h), preserve_range=True)
    else:
        image = skimage.transform.resize(image, (h * im_size / w, im_size), preserve_range=True)

    return image

def extract_aspect_box(image, aspect_ratio_box):
    """
    Extract a cropped region from the input image such that it matches a given aspect_ratio. The crop is
    extracted such that one dimension is preserved and the other is centrally cropped to fit the ratio.
    """
    image_h = image.shape[0]
    image_w = image.shape[1]

    aspect_ratio_image = get_aspect_ratio(image)

    if aspect_ratio_box >= aspect_ratio_image:
        box_h = round(image_w / aspect_ratio_box)
        image = image[image_h//2-box_h//2:image_h//2+box_h//2, ...]
    else:
        box_w = aspect_ratio_box * image_h
        image = image[:, image_w//2-box_w//2:image_w//2+box_w//2, ...]

    return image

def resize_image(image, size):
    """
    Just a wrapper for the skimage resize function
    """
    return skimage.transform.resize(image, size, preserve_range=True)

def get_aspect_ratio(image):
    """
    Returns the aspect ratio of an image
    """
    h, w, _ = image.shape
    return float(w) / h

def prepare_image(image, mean_values):
    """
    Get an image into the right format to feed the neural network.
    """

    # Shuffle axes to c01
    image = np.swapaxes(np.swapaxes(image, 1, 2), 0, 1)

    # Convert RGB to BGR
    image = image[::-1, :, :]

    # substract train sample mean
    image = image - mean_values

    return floatX(image[np.newaxis])

def deprocess(x, mean_values):
    """
    Convert the network output back into an image
    """
    x = np.copy(x[0])
    x += mean_values

    x = x[::-1]
    x = np.swapaxes(np.swapaxes(x, 0, 1), 1, 2)

    x = np.clip(x, 0, 255).astype('uint8')
    return x

max_stars_repo_name: baumgach/artistic-style-transfer

max_issues_repo_name: baumgach/artistic-style-transfer

max_forks_repo_name: baumgach/artistic-style-transfer


text:
"""Compute a Vertical Transect from an Extracted Polyline

usage:

    > python -O VerticalTransect.py

Necessary modules
-----------------
tputils
    Generic PyTecplot Utilities
numpy
    A general-purpose array-processing package
    https://pypi.org/project/numpy/

Description
-----------
To run this script we must first extract a Polyline geometry to be used in the transect.
Enable PyTecplot Connections via the Scripting menu to connect
to the running instance of Tecplot 360.

This script will connect to Tecplot 360, find the zone called Extracted Points
and use the XY points to define a vertical surface zone through the volume,
which defines the shape of the transect.
The solution variables are then interpolated from the solution data onto
the transect zone at each time step.

"""
import time
import numpy as np
import tecplot as tp
from tecplot.constant import *
import tputils

assert tp.version_info >= (0, 14), 'This script requires PyTecplot 0.14 or later'

tp.session.connect()

def get_first_zone_in_strand(strand, dataset):
    """Searches the zones by strand. Returns the first zone in the strand."""
    for z in dataset.zones():
        if z.strand == strand:
            return z
    return None

def plot_vertical_transect(zone, xvar, yvar):
    """Creates a new frame to plot the new Transect in a 2D Cartesian plot"""
    frame = tp.active_page().add_frame()
    tputils.attach_dataset(frame, zone.dataset)
    plot = frame.plot(PlotType.Cartesian2D)
    plot.activate()
    fmap = plot.fieldmap(zone)
    for f in plot.fieldmaps():
        f.show = False
    fmap.show = True
    plot.show_mesh = False
    plot.show_shade = False
    # Need to assign a contour variable before we can turn on the contour layer
    plot.contour(0).variable_index = 5  # FVCOM data, the first 5 variables are spatial
    plot.show_contour = True
    plot.show_scatter = False
    plot.axes.x_axis.variable = xvar
    plot.axes.y_axis.variable = yvar
    plot.axes.axis_mode = AxisMode.Independent
    # Fit the data first to get the ranges reasonable, then fine tune
    plot.view.fit()
    x_minmax = tputils.fieldmap_minmax(fmap, xvar)
    y_minmax = tputils.fieldmap_minmax(fmap, yvar)
    plot.axes.x_axis.min = x_minmax[0]
    plot.axes.x_axis.max = x_minmax[1]
    plot.axes.y_axis.min = y_minmax[0]
    plot.axes.y_axis.max = y_minmax[1]
    tp.macro.execute_command('$!Linking BetweenFrames{LinkSolutionTime = Yes}')
    tp.macro.execute_command('''$!PropagateLinking
      LinkType = BetweenFrames
      FrameCollection = All''')
    tp.macro.execute_extended_command(command_processor_id='Multi Frame Manager',
                                      command='TILEFRAMESHORIZ')
    return plot

def make_transect_zone(x_vals, y_vals, source_strand, dataset):
    """Creates a vertical transect from Siglev values"""
    dist_var = dataset.variable("distance")
    if dist_var is None:
        dist_var = dataset.add_variable("distance")
    dd = np.sqrt((x_vals[1:] - x_vals[:-1])**2 + (y_vals[1:] - y_vals[:-1])**2)
    dist = np.cumsum(dd)
    dist = np.concatenate([[0], dist])

    siglev_vals = np.unique(get_first_zone_in_strand(
        source_strand, dataset).values("siglev").as_numpy_array())

    iMax = len(x_vals)
    jMax = len(siglev_vals)

    X, Z = np.meshgrid(x_vals, siglev_vals)
    Y, Z = np.meshgrid(y_vals, siglev_vals)
    D, Z = np.meshgrid(dist, siglev_vals)

    zone_name = "Transect"
    result = dataset.add_ordered_zone(zone_name, (iMax, jMax, 1), locations=[
                                      ValueLocation.Nodal] * dataset.num_variables)
    result.values("x")[:] = X
    result.values("y")[:] = Y
    result.values("distance")[:] = D
    result.values("siglev")[:] = Z
    return result

def get_xy_points(dataset):
    """Looks for a zone called Extracted Points and imports the XY values into an array
    Could be edited to use a CSV or text file to define the Transect points."""
    transect_points_zone = dataset.zone("Extracted Points")
    x_vals = transect_points_zone.values("X").as_numpy_array()
    y_vals = transect_points_zone.values("Y").as_numpy_array()
    return x_vals, y_vals

def main():
    frame = tp.active_frame()
    threed_plot = tp.active_frame().plot(PlotType.Cartesian3D)
    dataset = tp.active_frame().dataset

    # FVCOM volume data should be strand #1
    source_strand = 1
    source_zones = [z for z in dataset.zones() if z.strand == source_strand]
    solution_time_to_zone = dict()
    for z in source_zones:
        solution_time_to_zone[z.solution_time] = z

    # This will interpolate nearly all the variables in the dataset.  If you're
    # only interested in a couple variables, then you can modify this list. This will
    # speed up the process, but could lead to confusing results since the vertical profile
    # will have "uninterpolated" values for any unspecified variables.  Buyer beware.
    variables_to_interpolate = list(dataset.variables())
    variables_to_interpolate.remove(dataset.variable("x"))
    variables_to_interpolate.remove(dataset.variable("y"))
    variables_to_interpolate.remove(dataset.variable("siglev"))
    # Might not be in the dataset yet
    dist_var = dataset.variable("distance")
    if dist_var is not None:
        variables_to_interpolate.remove(dist_var)

    begin = time.time()
    transect_base_zone = None
    orig_z_axis = threed_plot.axes.z_axis.variable

    # Save off the current view of the 3D frame and then set the Z-axis to siglev
    # so we can to a proper linear interpolation
    tp.macro.execute_command("$!View Copy")
    threed_plot.axes.z_axis.variable = dataset.variable("siglev")

    # For each time step we create an IJ ordered zone with XY points defined by x_vals & y_vals. The
    # z-values are the unique 'siglev' values.  We then interpolate from the volume data to this
    # IJ-ordered zone to create the transect.
    #
    # You must have enough XY points to capture the grid resolution of the volume dataset.
    x_vals, y_vals = get_xy_points(dataset)
    new_strand = tputils.max_strand(dataset) + 1

    for t in dataset.solution_times:
        start = time.time()
        if not transect_base_zone:
            transect_base_zone = make_transect_zone(x_vals, y_vals, source_strand, dataset)
            transect_base_zone.strand = new_strand
            transect_zone = transect_base_zone
        else:
            transect_zone = dataset.copy_zones([transect_base_zone])[0]
        transect_zone.solution_time = t
        source_zone = solution_time_to_zone[t]
        tp.data.operate.interpolate_linear(
            transect_zone, source_zones=source_zone, variables=variables_to_interpolate)
        print("Transect at t=%f complete: " % (t), time.time() - start)

    # Restore the 3D plot to its original view
    threed_plot.axes.z_axis.variable = orig_z_axis

    # Plot the results
    transect_plot = plot_vertical_transect(
        transect_base_zone, dataset.variable("distance"), dataset.variable("z"))

    # Set the contour color to match the 3D plot
    cont_group = threed_plot.fieldmap(source_zones[0]).contour.flood_contour_group_index
    transect_plot.contour(cont_group).variable_index = threed_plot.contour(
        cont_group).variable_index
    transect_plot.contour(cont_group).colormap_name = threed_plot.contour(cont_group).colormap_name
    transect_plot.fieldmap(transect_base_zone).contour.flood_contour_group_index = cont_group
    transect_plot.contour(cont_group).levels.reset_levels(*threed_plot.contour(cont_group).levels)

    frame.activate()
    tp.macro.execute_command("$!View Paste")

    print("Total time: ", time.time() - begin)

with tp.session.suspend():
    main()


text:
# Manual Linear Regression +optimization in TF

import numpy as np
import tensorflow as tf

def generate_dataset():
    x_batch = np.linspace(0, 2, 100)
    y_batch = 1.5 * x_batch + 0.5
    return x_batch, y_batch

def mean_squared_error( Y , y_pred ):
    return tf.reduce_mean( tf.square( y_pred - Y ) )

def mean_squared_error_deriv( Y , y_pred ):
    return tf.reshape( tf.reduce_mean( 2 * ( y_pred - Y ) ) , [ 1 , 1 ] )

def h ( X , weights , bias ):
    return tf.multiply( X , weights) + bias

x, y = generate_dataset()

num_features = 1
weight = tf.random.normal( ( num_features , 1 ) )
bias = tf.random.normal( ( num_features , 1 ) )

epochs_plot = list()
loss_plot = list()

num_epochs = 100
learning_rate = 0.01

for i in range( num_epochs ) :

    epoch_loss = list()
    for b in range( x.shape[0] ):

        output = h( x[b] , weight , bias )

        # Calculate loss and gradients
        loss = epoch_loss.append( mean_squared_error( y[b] , output ).numpy() )
        dJ_dH = mean_squared_error_deriv( y[b] , output)
        dH_dW = x[b]
        dJ_dW = tf.reduce_mean( dJ_dH * dH_dW )
        dJ_dB = tf.reduce_mean( dJ_dH )

        # Manual Optimization step
        weight -= ( learning_rate * dJ_dW )
        bias -= ( learning_rate * dJ_dB )

    loss = np.array( epoch_loss ).mean()
    epochs_plot.append( i + 1 )
    loss_plot.append( loss )

    print( 'Loss is {}'.format( loss ) )

print(weight, bias)

max_stars_repo_name: sunjaeyoon/LearningPytorch

max_issues_repo_name: sunjaeyoon/LearningPytorch

max_forks_repo_name: sunjaeyoon/LearningPytorch


text:
import logging
import os
import random
import shutil
import numpy as np

import sys
from pathlib import Path

import h5py

sys.path.append('~/carla/CARLA_0.9.6/carla/dist/carla-0.9.6-py3.5-linux-x86_64.egg')

import imageio

import image_converter
from config import SimConfig
from recordings import Recording
from sides import Side, SideMap

from carla_constants import *

import carla

class CarlaSim:
    def __init__(self,
                 config: SimConfig,
                 base_output_folder: str = '/data/carla/',
                 overwrite: bool = False,
                 seed=123,
                 ticks_per_frame=4,
                 car_idx=None,
                 host='localhost', port='2000'):

        random.seed(seed)

        self.config = config
        self.ticks_per_frame = ticks_per_frame
        self.num_peds = config.pedestrians
        self.num_cars = config.cars
        self.host = host
        self.port = port
        self.rgb_cameras = SideMap()
        self.depth_cameras = SideMap()
        self.ticks = 0

        self.output_folder = base_output_folder + self.config.folder_name + '/'

        self.frame_output_folder = self.output_folder + "raw/frames/"

        self.recording = Recording(base_output_folder, config)

        if os.path.exists(self.output_folder):
            if overwrite:
                shutil.rmtree(self.output_folder)
            else:
                print(
                    f"Output folder {self.output_folder} already exists and "
                    f"'overwrite' is false.")
                raise FileExistsError
        else:
            Path(self.output_folder).mkdir(exist_ok=True, parents=True)

            Path(self.frame_output_folder).mkdir(exist_ok=True, parents=True)
        with open(self.output_folder + "raw/seed.txt", 'x') as seed_file:
            seed_file.write(str(seed))

        self.pose_file = h5py.File(self.output_folder + "raw/pose.hdf5", 'w')
        self.abs_pose = self.pose_file.create_dataset("abs_pose", (0, 6), 'float32', maxshape=(None, 6))
        self.rel_pose = self.pose_file.create_dataset("rel_pose", (0, 6), 'float32', maxshape=(None, 6))
        self.start_rel_pose = self.pose_file.create_dataset("start_rel_pose", (0, 6), 'float32', maxshape=(None, 6))
        self.last_pose = None
        self.start_pose = None

        self.frames = 0

        try:
            # First of all, we need to create the client that will send the
            # requests
            # to the simulator. Here we'll assume the simulator is accepting
            # requests in the localhost at port 2000.
            self.client = carla.Client(self.host, int(self.port))
            self.client.set_timeout(2.0)

            # Once we have a client we can retrieve the world that is currently
            # running.
            # print(self.client.get_available_maps())
            self.world = self.client.load_world(self.config.city.map)

            # The world contains the list blueprints that we can use for
            # adding new
            # actors into the simulation.

            print(f"Spawning {self.num_cars} cars")

            blueprints = self.world.get_blueprint_library().filter("vehicle.*")
            blueprintsWalkers = self.world.get_blueprint_library().filter(
                "walker.pedestrian.*")

            blueprints = [x for x in blueprints if
                          int(x.get_attribute('number_of_wheels')) == 4]
            blueprints = [x for x in blueprints if
                          not x.id.endswith('isetta')]
            blueprints = [x for x in blueprints if
                          not x.id.endswith('carlacola')]

            self.cars = []

            spawn_points = self.world.get_map().get_spawn_points()
            number_of_spawn_points = len(spawn_points)

            if self.num_cars <= number_of_spawn_points:
                random.shuffle(spawn_points)
            else:
                raise ValueError(
                    f"Number of cars {self.num_cars} greater than number of "
                    f"spawn "
                    f"points {number_of_spawn_points}")

            SpawnActor = carla.command.SpawnActor
            SetAutopilot = carla.command.SetAutopilot
            FutureActor = carla.command.FutureActor

            batch = []
            for n, transform in enumerate(spawn_points):

                if n >= self.num_cars:
                    break

                blueprint = random.choice(blueprints)
                if blueprint.has_attribute('color'):
                    color = random.choice(
                        blueprint.get_attribute('color').recommended_values)
                    blueprint.set_attribute('color', color)
                if blueprint.has_attribute('driver_id'):
                    driver_id = random.choice(
                        blueprint.get_attribute('driver_id').recommended_values)
                    blueprint.set_attribute('driver_id', driver_id)
                blueprint.set_attribute('role_name', 'autopilot')
                batch.append(SpawnActor(blueprint, transform).then(
                    SetAutopilot(FutureActor, True)))

            for response in self.client.apply_batch_sync(batch):
                if response.error:
                    logging.error(response.error)
                else:
                    # print("spawned car", response.actor_id)
                    self.cars.append(response.actor_id)

            self.cars = [self.world.get_actor(x) for x in self.cars]

            # -----------  WALKERS  -----------

            print(f"Spawning {self.num_peds} walkers")

            spawn_points = []
            for i in range(self.num_peds):
                spawn_point = carla.Transform()
                loc = self.world.get_random_location_from_navigation()
                if loc is not None:
                    spawn_point.location = loc
                    spawn_points.append(spawn_point)

            walkers_list = []
            all_id = []

            # Spawn walkers
            batch = []
            for spawn_point in spawn_points:
                walker_bp = random.choice(blueprintsWalkers)
                # set as not invencible
                if walker_bp.has_attribute('is_invincible'):
                    walker_bp.set_attribute('is_invincible', 'false')
                batch.append(SpawnActor(walker_bp, spawn_point))
            results = self.client.apply_batch_sync(batch, True)
            for i in range(len(results)):
                if results[i].error:
                    logging.error(results[i].error)
                else:
                    # print("spawned walker", results[i].actor_id)
                    walkers_list.append({"id": results[i].actor_id})

            batch = []
            walker_controller_bp = self.world.get_blueprint_library().find(
                'controller.ai.walker')
            for i in range(len(walkers_list)):
                batch.append(SpawnActor(walker_controller_bp, carla.Transform(),
                                        walkers_list[i]["id"]))
            results = self.client.apply_batch_sync(batch, True)
            for i in range(len(results)):
                if results[i].error:
                    logging.error(results[i].error)
                else:
                    walkers_list[i]["con"] = results[i].actor_id

            # 4. we put altogether the walkers and controllers id to get the
            # objects from their id
            for i in range(len(walkers_list)):
                all_id.append(walkers_list[i]["con"])
                all_id.append(walkers_list[i]["id"])
            all_actors = self.world.get_actors(all_id)

            # wait for a tick to ensure client receives the last transform of
            # the walkers we have just created
            self.world.wait_for_tick()

            # 5. initialize each controller and set target to walk to (list
            # is [controler, actor, controller, actor ...])
            for i in range(0, len(all_id), 2):
                # start walker
                all_actors[i].start()
                # set walk to random point
                all_actors[i].go_to_location(
                    self.world.get_random_location_from_navigation())
                # random max speed
                all_actors[i].set_max_speed(
                    1 + random.random())  # max speed between 1 and 2 (
                # default is 1.4 m/s)

            self.pedestrians = [x['id'] for x in walkers_list]
            self.pedestrians = [self.world.get_actor(x) for x in
                                self.pedestrians]

            blueprint_library = self.world.get_blueprint_library()

            vehicles = blueprint_library.filter('vehicle')

            vehicles = [v for v in vehicles if v.get_attribute('number_of_wheels').as_int() >= 4]

            # Add main car
            if car_idx is None:
                bp = random.choice(vehicles)
            else:
                bp = vehicles[car_idx]

            if bp.has_attribute('color'):
                color = random.choice(
                    bp.get_attribute('color').recommended_values)
                bp.set_attribute('color', color)

            transform = random.choice(self.world.get_map().get_spawn_points())

            # So let's tell the world to spawn the vehicle.
            vehicle = self.world.spawn_actor(bp, transform)

            print('created %s' % vehicle.type_id)

            vehicle.set_autopilot(True)

            self.car = vehicle
            self.world.get_spectator().set_transform(self.car.get_transform())

            self.add_sensors()
            self.world.set_weather(config.weather.params)

            settings = self.world.get_settings()
            settings.synchronous_mode = True
            settings.fixed_delta_seconds = 0.05
            settings.no_rendering_mode = False
            self.world.apply_settings(settings)

        except Exception as e:
            print("Error starting server")
            raise e

    def add_sensors(self):
        blueprint_library = self.world.get_blueprint_library()

        for side in list(Side):
            camera_bp = blueprint_library.find('sensor.camera.rgb')
            camera_bp.set_attribute('image_size_x', str(IMAGE_WIDTH))
            camera_bp.set_attribute('image_size_y', str(IMAGE_HEIGHT))
            camera_bp.set_attribute('fov', str(FOV))

            camera_transform = carla.Transform(carla.Location(x=1.5, z=2.4),
                                               side.rotation())
            camera = self.world.spawn_actor(camera_bp, camera_transform,
                                            attach_to=self.car)
            # print('created %s' % camera.type_id)
            self.rgb_cameras[side] = camera

        for side in list(Side):
            camera_bp = blueprint_library.find('sensor.camera.depth')
            camera_bp.set_attribute('image_size_x', str(IMAGE_WIDTH))
            camera_bp.set_attribute('image_size_y', str(IMAGE_HEIGHT))
            camera_bp.set_attribute('fov', str(FOV))

            camera_transform = carla.Transform(carla.Location(x=1.5, z=2.4),
                                               side.rotation())
            camera = self.world.spawn_actor(camera_bp, camera_transform,
                                            attach_to=self.car)
            # print('created %s' % camera.type_id)
            self.depth_cameras[side] = camera

    def tick(self):
        frame = self.world.tick()

        car_loc = self.car.get_transform()
        vec = car_loc.get_forward_vector()
        loc = car_loc.location

        new_loc = carla.Location(loc.x - vec.x * 12, loc.y - vec.y * 12,
                                 loc.z + 3)

        self.world.get_spectator().set_transform(
            carla.Transform(new_loc, car_loc.rotation))

        do_save = self.ticks % self.ticks_per_frame == 0
        if do_save:
            images = self.rgb_cameras.pop(frame)
            for k, v in images.items():
                v.save_to_disk(
                    self.frame_output_folder + f"{k.name.lower()}_"
                                               f"{int(self.frames)}_rgb.png")
            images = self.depth_cameras.pop(frame)
            for k, v in images.items():
                norm_depth = image_converter.depth_to_array(v)

                imageio.imwrite(
                    self.frame_output_folder + f"{k.name.lower()}_"
                                               f"{int(self.frames)}_depth.png",
                    (norm_depth * DEPTH_MULTIPLIER).astype('uint16'))

            orientation = np.array([vec.x, vec.y, vec.z])
            pose = np.array([loc.x, loc.y, loc.z])
            pose = np.concatenate([pose, orientation])

            self.abs_pose.resize(self.frames + 1, axis=0)
            self.rel_pose.resize(self.frames + 1, axis=0)
            self.start_rel_pose.resize(self.frames + 1, axis=0)

            self.abs_pose[-1] = pose

            if self.last_pose is None:
                self.rel_pose[-1] = np.zeros((6,), 'float32')
                self.start_rel_pose[-1] = np.zeros((6,), 'float32')

                self.start_pose = pose
            else:
                self.rel_pose[-1] = pose - self.last_pose
                self.start_rel_pose[-1] = pose - self.start_pose

            self.last_pose = pose

            self.frames += 1

        self.ticks += 1

        return do_save

    def end(self):
        for a in self.cars:
            a.destroy()
        for a in self.pedestrians:
            a.destroy()

        for a in self.rgb_cameras.sides.values():
            a[0].destroy()
        for a in self.depth_cameras.sides.values():
            a[0].destroy()

        self.car.destroy()

max_stars_repo_name: rnett/CARLASim

max_issues_repo_name: rnett/CARLASim

max_forks_repo_name: rnett/CARLASim


text:
# Copyright 2020 JD.com, Inc. Galileo Authors. 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 os
import pytest
import tensorflow as tf
from galileo.tests import expected_data
from galileo.tests.utils import numpy_equal_unique, numpy_equal
from galileo.tf import ops

valid_ids = ('global', 'valid_zero', 'valid_one', 'valid_multi')
valid_topk_ids = ('global', 'valid_one', 'valid_multi')
invalid_ids = ('no_neighbor', 'invalid_types')

valid_params = (
    ([1001], [], [1000, 1003, 1004, 1005], [3.5, 5.5]),
    ([1001], [0], [1003, 1004, 1005], [3.5]),
    ([1001], [1], [1000], [5.5]),
    ([1000, 1002, 1006, 1009], [0, 1], [1001, 1001, 1000, 1000], [3.5, 5.5]),
)

valid_topk_params = (
    ([1001], [], [1000], [5.5]),
    ([1001], [1], [1000], [5.5]),
    ([1000, 1002, 1006, 1009], [0, 1], [1001, 1001, 1000,
                                        1000], [3.5, 3.5, 5.5, 5.5]),
)

invalid_params = (
    pytest.param([1003, 1004, 1005], []),
    pytest.param([1001], [-1, 5]),
)

@pytest.mark.parametrize('vertex,types,expected_nbr,expected_weight',
                         valid_params,
                         ids=valid_ids)
def test_valid_sample_neighbor_with_weight(prepare_tf_env, vertex, types,
                                           expected_nbr, expected_weight):
    count = 100
    is_weight = True
    neighbors, weights = ops.sample_neighbors(vertex, types, count, is_weight)
    assert numpy_equal([len(vertex), count], tf.shape(neighbors).numpy())
    assert numpy_equal([len(vertex), count], tf.shape(weights).numpy())
    assert numpy_equal_unique(expected_nbr, neighbors.numpy())
    assert numpy_equal_unique(expected_weight, weights.numpy())

@pytest.mark.parametrize('vertex,types', invalid_params, ids=invalid_ids)
def test_invalid_sample_neighbor_with_weight(prepare_tf_env, vertex, types):
    count = 100
    is_weight = True
    with pytest.raises(Exception):
        ops.sample_neighbors(vertex, types, count, is_weight)

@pytest.mark.parametrize('vertex,types,expected_nbr,expected_weight',
                         valid_params,
                         ids=valid_ids)
def test_valid_sample_neighbor_without_weight(prepare_tf_env, vertex, types,
                                              expected_nbr, expected_weight):
    vertex_tensor = tf.constant(vertex, dtype=tf.int64)
    types_tensor = tf.constant(types, dtype=tf.uint8)
    count = 100
    is_weight = False
    neighbors = ops.sample_neighbors(vertex_tensor, types_tensor, count,
                                     is_weight)[0]
    assert numpy_equal([len(vertex), count], tf.shape(neighbors).numpy())
    assert numpy_equal_unique(expected_nbr, neighbors.numpy())

@pytest.mark.parametrize('vertex,types', invalid_params, ids=invalid_ids)
def test_invalid_sample_neighbor_without_weight(prepare_tf_env, vertex, types):
    vertex_tensor = tf.constant(vertex, dtype=tf.int64)
    types_tensor = tf.constant(types, dtype=tf.uint8)
    count = 100
    is_weight = False
    with pytest.raises(Exception):
        ops.sample_neighbors(vertex_tensor, types_tensor, count, is_weight)[0]

@pytest.mark.parametrize('vertex,types,expected_nbr,expected_weight',
                         valid_topk_params,
                         ids=valid_topk_ids)
def test_valid_topk_neighbor_with_weight(prepare_tf_env, vertex, types,
                                         expected_nbr, expected_weight):
    count = 1
    is_weight = True
    neighbors = ops.get_topk_neighbors(vertex, types, count, is_weight)
    assert numpy_equal([len(vertex), 1], tf.shape(neighbors[0]).numpy())
    assert numpy_equal([len(vertex), 1], tf.shape(neighbors[1]).numpy())
    assert numpy_equal(expected_nbr, tf.reshape(neighbors[0], [-1]).numpy())
    assert numpy_equal(expected_weight, tf.reshape(neighbors[1], [-1]).numpy())

@pytest.mark.parametrize('vertex,types', invalid_params, ids=invalid_ids)
def test_invalid_topk_neighbor_with_weight(prepare_tf_env, vertex, types):
    count = 2
    is_weight = True
    with pytest.raises(Exception):
        ops.get_topk_neighbors(vertex, types, count, is_weight)

@pytest.mark.parametrize('vertex,types,expected_nbr,expected_weight',
                         valid_topk_params,
                         ids=valid_topk_ids)
def test_valid_topk_neighbor_without_weight(prepare_tf_env, vertex, types,
                                            expected_nbr, expected_weight):
    vertex_tensor = tf.constant(vertex, dtype=tf.int64)
    types_tensor = tf.constant(types, dtype=tf.uint8)
    count = 1
    is_weight = False
    neighbors = ops.get_topk_neighbors(vertex_tensor, types_tensor, count,
                                       is_weight)
    assert numpy_equal([len(vertex), 1], tf.shape(neighbors[0]).numpy())
    assert numpy_equal(expected_nbr, tf.reshape(neighbors[0], [-1]).numpy())

@pytest.mark.parametrize('vertex,types', invalid_params, ids=invalid_ids)
def test_invalid_topk_neighbor_without_weight(prepare_tf_env, vertex, types):
    count = 2
    is_weight = False
    vertex_tensor = tf.constant(vertex, dtype=tf.int64)
    types_tensor = tf.constant(types, dtype=tf.uint8)
    with pytest.raises(Exception):
        ops.get_topk_neighbors(vertex_tensor, types_tensor, count, is_weight)

@pytest.mark.parametrize('vertex,types,expected_nbr,expected_weight',
                         valid_params,
                         ids=valid_ids)
def test_valid_full_neighbor_with_weight(prepare_tf_env, vertex, types,
                                         expected_nbr, expected_weight):
    is_weight = True
    neighbors = ops.get_full_neighbors(vertex, types, is_weight)
    assert 3 == len(neighbors)
    assert numpy_equal_unique(expected_nbr, neighbors[0].numpy())
    assert numpy_equal_unique(expected_weight, neighbors[1].numpy())

@pytest.mark.parametrize('vertex,types', invalid_params, ids=invalid_ids)
def test_invalid_full_neighbor_with_weight(prepare_tf_env, vertex, types):
    vertex_tensor = tf.constant(vertex, dtype=tf.int64)
    types_tensor = tf.constant(types, dtype=tf.uint8)
    is_weight = True
    with pytest.raises(Exception):
        ops.get_full_neighbors(vertex_tensor, types_tensor, is_weight)

@pytest.mark.parametrize('vertex,types,expected_nbr,expected_weight',
                         valid_params,
                         ids=valid_ids)
def test_valid_full_neighbor_without_weight(prepare_tf_env, vertex, types,
                                            expected_nbr, expected_weight):
    vertex_tensor = tf.constant(vertex, dtype=tf.int64)
    types_tensor = tf.constant(types, dtype=tf.uint8)
    is_weight = False
    neighbors = ops.get_full_neighbors(vertex_tensor, types_tensor, is_weight)
    assert 2 == len(neighbors)
    assert numpy_equal_unique(expected_nbr, neighbors[0].numpy())

@pytest.mark.parametrize('vertex,types', invalid_params, ids=invalid_ids)
def test_invalid_full_neighbor_without_weight(prepare_tf_env, vertex, types):
    vertex_tensor = tf.constant(vertex, dtype=tf.int64)
    types_tensor = tf.constant(types, dtype=tf.uint8)
    is_weight = False
    with pytest.raises(Exception):
        ops.get_full_neighbors(vertex_tensor, types_tensor, is_weight)

def test_empty_inputs(prepare_tf_env):
    vertex_tensor = tf.constant(1, dtype=tf.int64, shape=(0, ))
    types_tensor = tf.constant([0], dtype=tf.uint8)

    res = ops.sample_neighbors(vertex_tensor, types_tensor, 5, False)
    assert 1 == len(res)
    assert res[0].dtype == tf.int64
    assert res[0].shape[0] == 0
    assert res[0].shape[1] == 5

    res = ops.get_topk_neighbors(vertex_tensor, types_tensor, 2, False)
    assert 1 == len(res)
    assert res[0].dtype == tf.int64
    assert res[0].shape[0] == 0
    assert res[0].shape[1] == 2

    res = ops.get_full_neighbors(vertex_tensor, types_tensor, False)
    assert 2 == len(res)
    assert res[0].dtype == tf.int64
    assert res[0].shape[0] == 0
    assert res[1].dtype == tf.int32
    assert res[1].shape[0] == 0
    assert res[1].shape[1] == 2

    res = ops.sample_neighbors(vertex_tensor, types_tensor, 5, True)
    assert 2 == len(res)
    assert res[0].dtype == tf.int64
    assert res[0].shape[0] == 0
    assert res[0].shape[1] == 5
    assert res[1].dtype == tf.float32
    assert res[1].shape[0] == 0
    assert res[1].shape[1] == 5

    res = ops.get_topk_neighbors(vertex_tensor, types_tensor, 2, True)
    assert 2 == len(res)
    assert res[0].dtype == tf.int64
    assert res[0].shape[0] == 0
    assert res[0].shape[1] == 2
    assert res[1].dtype == tf.float32
    assert res[1].shape[0] == 0
    assert res[1].shape[1] == 2

    res = ops.get_full_neighbors(vertex_tensor, types_tensor, True)
    assert 3 == len(res)
    assert res[0].dtype == tf.int64
    assert res[0].shape[0] == 0
    assert res[1].dtype == tf.float32
    assert res[1].shape[0] == 0
    assert res[2].dtype == tf.int32
    assert res[2].shape[0] == 0
    assert res[2].shape[1] == 2

max_issues_repo_name: Hacky-DH/galileo

max_forks_repo_name: Hacky-DH/galileo


text:
from qcor import *

# Simple 1-qubit demonstration of the Quatum Imaginary Time Evolution/QLanczos algorithm
# Reference: https://www.nature.com/articles/s41567-019-0704-4
# Target H = 1/sqrt(2)(X + Z)
# Expected minimum value: -1.0
H = 0.7071067811865475 * X(0) + 0.7071067811865475 * Z(0)

# See Fig. 2 (e) of https://www.nature.com/articles/s41567-019-0704-4
# Horizontal axis: 0 -> 2.5
# The number of Trotter steps
nbSteps = 25

# The Trotter step size
stepSize = 0.1

# Create the problem model
problemModel = QuaSiMo.ModelFactory.createModel(H)

# Create the QITE workflow
workflow = QuaSiMo.getWorkflow('qite', {'steps': nbSteps, 'step-size': stepSize})

# Execute and print the result
result = workflow.execute(problemModel)

# Final energy:
energy = result['energy']
print('Ground state energy =', energy)

# Plot energy value along QITE steps
td_energy_vals = result['exp-vals']
import matplotlib.pyplot as plt
import numpy as np
plt.plot(np.arange(0, nbSteps + 1) * stepSize, td_energy_vals, 'ro-', label = 'QITE')
plt.legend()
plt.grid()
plt.savefig('qite.pdf')

max_stars_repo_name: ausbin/qcor

max_issues_repo_name: ausbin/qcor

max_forks_repo_name: ausbin/qcor


text:
from datetime import date, timedelta
from datetime import datetime as dt
from matplotlib.font_manager import findfont, FontProperties
from PIL import Image
from PIL import ImageFont
from PIL import ImageDraw
import numpy as np
import threading
import pickle
import time
import json
import math
import sys
import os
import matplotlib
matplotlib.use('Agg')
import matplotlib.pyplot as plt

PKL1 = 'data/real_data_1.pkl'
PKL2 = 'data/real_data_2.pkl'

backends = ['ibmqx4', 'ibmqx5']

def get_token(path):
    with open(path) as jsn:
        data = json.load(jsn)
    return data['token']

def load_data():
    data = None
    try:
        with open(PKL1, 'rb') as f:
            data = pickle.load(f)
    except (pickle.UnpicklingError, EOFError) as e:
        with open(PKL2, 'rb') as f:
            data = pickle.load(f)
    return data

def plot_pending_jobs(backend):
    data = load_data()

    times = sorted([x[0] for x in data])
    pending_jobs = [[y for y in x[1] if y['backend'] == backend][0]['pending_jobs']
                    for x in sorted(data, key=lambda x: x[0])]

    plt.figure(figsize=(11, 5))
    plt.grid(True, zorder=5)
    plt.fill_between(times, pending_jobs, color='brown')

    # New xticks.
    locs, labels = plt.xticks()
    new_ticks = [dt.fromtimestamp(x).strftime('%H:%M') for x in locs]
    plt.xticks(locs[1:-1], new_ticks[1:-1], rotation=0, fontsize=15)
    plt.yticks(fontsize=15)

    plt.ylim(0, math.ceil(max(pending_jobs)) + 1)

    plt.title('Local time of bot: {}'.format(dt.fromtimestamp(time.time()).strftime('%Y, %b %d, %H:%M')),
              fontsize=15)
    plt.xlabel('Time', fontsize=15)
    plt.ylabel('# of pending jobs', fontsize=15)

    filename = 'tmp/{}.png'.format(backend)
    plt.savefig(filename, bbox_inches='tight')
    plt.close()

def plot_calibration(backend, api):
    full_info = api.backend_calibration(backend=backend)

    # Info.
    N_qubits = len(full_info['qubits'])
    qubits = [full_info['qubits'][qub]['name'] for qub in range(N_qubits)]
    readout_error = [full_info['qubits'][qub]['readoutError']['value'] for qub in range(N_qubits)]
    readout_error = np.array([readout_error])

    ###################
    # Multiqubit error
    multi_qubit_gates = [full_info['multiQubitGates'][qub]['qubits'] for qub in range(N_qubits)]
    multi_qubit_error = [full_info['multiQubitGates'][qub]['gateError']['value'] for qub in range(N_qubits)]

    # creating gate error matrix.
    error_matrix = np.zeros((N_qubits, N_qubits))

    error_matrix[:, :] = None
    for i in range(len(multi_qubit_gates)):
        gate = multi_qubit_gates[i]
        qub1, qub2 = gate[0], gate[1]
        error_matrix[qub1][qub2] = multi_qubit_error[i]
        error_matrix[qub2][qub1] = multi_qubit_error[i]
        error_matrix[i][i] = readout_error[0][i]

    plt.close()
    plt.figure(figsize=(6, 6))
    plt.matshow(error_matrix, cmap='Reds', fignum=1)
    # plt.title('Backend: {}, Two qubit gate errors,\n last calibration: {}'.format(backend, last_update), fontsize=15)
    plt.title('Two qubit gate errors', fontsize=15)
    # Placing actual values in the matshow plot
    for (i, j), value in np.ndenumerate(error_matrix):
        if not np.isnan(value):
            plt.text(j, i, '{:0.2f}'.format(value), ha='center', va='center')
    # Formatting axes
    plt.yticks(np.arange(N_qubits), qubits)
    plt.xticks(np.arange(N_qubits), qubits)
    plt.autoscale(axis='both', tight=True)
    plt.savefig('tmp/{}_multiqubut_err.png'.format(backend), bbox_inches='tight')
    plt.close()

def create_statistics(backend, api):
    plot_calibration(backend, api)
    plot_pending_jobs(backend)

    img_merror = Image.open('tmp/{}_multiqubut_err.png'.format(backend), 'r')
    img_merror_w, img_merror_h = img_merror.size

    img_jobs = Image.open('tmp/{}.png'.format(backend), 'r')
    img_jobs_w, img_jobs_h = img_jobs.size

    img_background = Image.new('RGBA', (img_merror_w + img_jobs_w, img_merror_h + 10),
                               (255, 255, 255, 255))
    img_background_w, img_background_h = img_background.size

    #############
    # Load Logos.
    img_qiskit = Image.open('res/qiskit-logo.png', 'r')
    factor = 5
    img_qiskit = img_qiskit.resize((img_qiskit.size[0] // factor, img_qiskit.size[1] // factor))
    img_qiskit_w, img_qiskit_h = img_qiskit.size

    img_rqc = Image.open('res/rqc.jpg', 'r')
    factor = 10
    img_rqc = img_rqc.resize((img_rqc.size[0] // factor, img_rqc.size[1] // factor))
    img_rqc_w, img_rqc_h = img_rqc.size
    ####

    ##############
    # Paste Plots.
    displacement_h = 10
    offset = (0, displacement_h)
    img_background.paste(img_merror, offset)

    displacement_h = 50
    offset = (img_merror_w, displacement_h)
    img_background.paste(img_jobs, offset)
    ####

    ##############
    # Paste Logos.
    displacement_h = 8
    displacement_w = 10
    offset = (img_background_w - img_qiskit_w - displacement_w, displacement_h - 6)
    img_background.paste(img_qiskit, offset)

    offset = (img_background_w - img_rqc_w - img_qiskit_w - displacement_w - 7, displacement_h)
    img_background.paste(img_rqc, offset)
    ###

    font_path = findfont(FontProperties(family=['sans-serif'], weight='bold'))

    full_info = api.backend_calibration(backend=backend)
    last_update = full_info['lastUpdateDate']
    last_update = dt.strptime(last_update, "%Y-%m-%dT%H:%M:%S.000Z").timestamp()
    last_update = dt.fromtimestamp(last_update).strftime('%Y, %b %d, %H:%M')
    draw = ImageDraw.Draw(img_background)
    font = ImageFont.truetype(font_path, 16)
    title = 'IBMQ Backend: {},\nLast calibration: {}'.format(backend, last_update)
    draw.text((img_merror_w + 45, 10), title, (0, 0, 0), font=font)

    img_background

    img_background.save('tmp/{}_to_send.png'.format(backend))

    img_background.close()
    img_rqc.close()
    img_merror.close()
    img_jobs.close()

if __name__ == '__main__':
    raise RuntimeError

max_stars_repo_name: akarazeevprojects/QuantumComputingBot

max_issues_repo_name: akarazeevprojects/QuantumComputingBot

max_forks_repo_name: akarazeevprojects/QuantumComputingBot