text:
############################################################
# Image augmentation utility functions file
############################################################

import cv2
import numpy as np
import numbers
import random
import collections

def get_2dshape(shape, *, zero=True):
    if not isinstance(shape, collections.Iterable):
        shape = int(shape)
        shape = (shape, shape)
    else:
        h, w = map(int, shape)
        shape = (h, w)
    if zero:
        minv = 0
    else:
        minv = 1

    assert min(shape) >= minv, 'invalid shape: {}'.format(shape)
    return shape

# function for cropping image at position generated with a specific size
def random_crop_pad_to_shape(img, crop_pos, crop_size, pad_label_value):
    h, w = img.shape[:2]
    start_crop_h, start_crop_w = crop_pos
    assert ((start_crop_h < h) and (start_crop_h >= 0))
    assert ((start_crop_w < w) and (start_crop_w >= 0))

    crop_size = get_2dshape(crop_size)
    crop_h, crop_w = crop_size

    img_crop = img[start_crop_h:start_crop_h + crop_h,
               start_crop_w:start_crop_w + crop_w, ...]

    img_, margin = pad_image_to_shape(img_crop, crop_size, cv2.BORDER_CONSTANT,
                                      pad_label_value)

    return img_, margin

# function to generate a random position for cropping
def generate_random_crop_pos(ori_size, crop_size):
    ori_size = get_2dshape(ori_size)
    h, w = ori_size

    crop_size = get_2dshape(crop_size)
    crop_h, crop_w = crop_size

    pos_h, pos_w = 0, 0

    if h > crop_h:
        pos_h = random.randint(0, h - crop_h + 1)

    if w > crop_w:
        pos_w = random.randint(0, w - crop_w + 1)

    return pos_h, pos_w

def pad_image_to_shape(img, shape, border_mode, value):
    margin = np.zeros(4, np.uint32)
    shape = get_2dshape(shape)
    pad_height = shape[0] - img.shape[0] if shape[0] - img.shape[0] > 0 else 0
    pad_width = shape[1] - img.shape[1] if shape[1] - img.shape[1] > 0 else 0

    margin[0] = pad_height // 2
    margin[1] = pad_height // 2 + pad_height % 2
    margin[2] = pad_width // 2
    margin[3] = pad_width // 2 + pad_width % 2

    img = cv2.copyMakeBorder(img, margin[0], margin[1], margin[2], margin[3],
                             border_mode, value=value)

    return img, margin

def random_scale(img, gt, scales):
    scale = random.choice(scales)
    sh = int(img.shape[0] * scale)
    sw = int(img.shape[1] * scale)
    img = cv2.resize(img, (sw, sh), interpolation=cv2.INTER_LINEAR)
    gt = cv2.resize(gt, (sw, sh), interpolation=cv2.INTER_NEAREST)

    return img, gt, scale

def random_mirror(img, gt):
    if random.random() >= 0.5:
        img = cv2.flip(img, 1)
        gt = cv2.flip(gt, 1)

    return img, gt

def normalize(img, mean, std):
    # pytorch pretrained model need the input range: 0-1
    img = img.astype(np.float32) / 255.0
    img = img - mean
    img = img / std

    return img

ext: py

reason: import numpy


text:
import numpy as np
import datetime
import sys
import os
import json
import pandas as pd
import matplotlib.pyplot as plt
from nilmtk.dataset import DataSet
from collections import Counter
from joblib import parallel_backend
from sklearn.model_selection import KFold
from sklearn.ensemble import IsolationForest
from sklearn.neighbors import LocalOutlierFactor
from sklearn.svm import OneClassSVM

## LOAD DATA FROM DATASET ##

def loadApplicationDataFromDataset(application, settings):
    data = DataSet(settings["h5File"])
    df = None

    window = settings["window"]
    if settings["window"] != None:
        data.set_window(window[0],window[1])
    try:
        elecMeter = data.buildings[application[0]].elec
        df = elecMeter[application[1]].power_series_all_data()
    except:
        print(str(e))
        return None

    if df is None:
      return None

    if df.isnull().values.any():
        df.fillna(method='ffill')

    if settings["filterData"] == "rollingMedian":
        df = df.rolling(window=10,
                        center=True,
                        min_periods=1).median()
    elif settings["filterData"] == "hampel":
        df = hampel_filter(df, 10)

    power_values = np.arange(int(max(df))+1).reshape(-1, 1)
    x = df.astype('int').values.reshape(-1, 1)

    return {"power_val": power_values, "x": x}

def loadAggregateDataFromDataset(application, settings):
    data = DataSet(settings["h5File"])
    df = None

    window = settings["window"]
    if settings["window"] != None:
        data.set_window(window[0],window[1])
    try:
        elecMeter = data.buildings[application[0]].elec
        df = elecMeter.mains().power_series_all_data()
    except Exception as e:
        print(str(e))
        return None

    if df is None:
      return None

    if df.isnull().values.any():
        df.fillna(method='ffill')

    if settings["filterData"] == "rollingMedian":
        df = df.rolling(window=10,
                        center=True,
                        min_periods=1).median()
    elif settings["filterData"] == "hampel":
        df = hampelFilter(df, 10)

    power_values = np.arange(int(max(df))+1).reshape(-1, 1)
    x = df.astype('int').values.reshape(-1, 1)

    return {"power_val": power_values, "x": x}

## HAMPEL FILTER ##

def hampelFilter(df, windowSize, nSig=3):
    # https://towardsdatascience.com/outlier-detection-with-hampel-filter-85ddf523c73d
    k = 1.4826 # Gaussian distribution scale factor
    dfFiltered = df.copy()

    MAD = lambda x: np.median(np.abs(x - np.median(x)))
    rollMedian = df.rolling(window=2*windowSize,
                            center=True).median()
    rollMAD = k * df.rolling(window=2*windowSize,
                            center=True).apply(MAD, raw=False)
    diff = np.abs(df - rollMedian)

    idx = list(np.argwhere(diff > (nSig * rollMAD)).flatten())
    dfFiltered[idx] = rollMedian[idx]

    return dfFiltered

## ISOLATION FORESTS ###

def calcIsoForest(x_train, x_test):
    try:
        clf = IsolationForest(random_state=0,
                                n_jobs=6).fit(x_train)
        outliers = clf.predict(x_test)
    except:
        return None
    return outliers

## Local Outlier Factor

def calcLOF(x_train, x_test):
    try:
        clf = LocalOutlierFactor(n_neighbors=2,
                                    n_jobs=6)
        outliers = clf.fit_predict(x_test)
    except:
        return None
    return outliers

def storeResults(settings, metrics):
    if os.path.exists(settings["fileName"]):
        fileMode = 'a'
    else:
        fileMode = 'w'

    with open(settings["fileName"], fileMode) as outfile:
        if settings["jsonIndent"]:
            json.dump([settings, metrics],
                        outfile,
                        indent=1)
        else:
            json.dump([settings, metrics], outfile)
        outfile.write("\n")

## METRICS ###

def calcMetrics(results, settings):
    metrics = {}

    if "lof" in settings["algo"]:
        metrics["lof"] = {}
        metrics["lof"]["absApplOut"] = np.count_nonzero(results["outLOFTest"] == -1)
        metrics["lof"]["relApplOut"] = np.count_nonzero(results["outLOFTest"] == -1)/len(results["outLOFTest"])
        lofPowerOutlier = (results["testData"][(np.where(results["outLOFTest"] == -1))])[:,0]
        keys_values = dict(Counter(lofPowerOutlier)).items()
        metrics["lof"]["applOutOccur"] = {str(key): str(value) for key, value in keys_values}
        metrics["lof"]["applOutlier"] = (np.where(results["outLOFTest"] == -1))[0].tolist()

        metrics["lof"]["absPowOut"] = np.count_nonzero(results["outLOFPower"] == -1)
        metrics["lof"]["relPowOut"] = np.count_nonzero(results["outLOFPower"] == -1)/len(results["outLOFPower"])
        metrics["lof"]["powerOutVal"] = (results["powerData"][(np.where(results["outLOFPower"] == -1))])[:,0].tolist()
        metrics["lof"]["powerOutlier"] = (np.where(results["outLOFPower"] == -1))[0].tolist()

    if "if" in settings["algo"]:
        metrics["if"]  = {}
        metrics["if"]["absApplOut"] = np.count_nonzero(results["outIsoForTest"] == -1)
        metrics["if"]["relApplOut"] = np.count_nonzero(results["outIsoForTest"] == -1)/len(results["outIsoForTest"])
        ifPowerOutlier = (results["testData"][(np.where(results["outIsoForTest"] == -1))])[:,0]
        keys_values = dict(Counter(ifPowerOutlier)).items()
        metrics["if"]["applOutOccur"] = {str(key): str(value) for key, value in keys_values}
        metrics["if"]["applOutlier"] = (np.where(results["outIsoForTest"] == -1))[0].tolist()

        metrics["if"]["absPowOut"] = np.count_nonzero(results["outIsoForPower"] == -1)
        metrics["if"]["relPowOut"] = np.count_nonzero(results["outIsoForPower"] == -1)/len(results["outIsoForPower"])
        metrics["if"]["powerOutVal"] = (results["powerData"][(np.where(results["outIsoForPower"] == -1))])[:,0].tolist()
        metrics["if"]["powerOutlier"] = (np.where(results["outIsoForPower"] == -1))[0].tolist()

    return metrics

def calcOutliers(data, settings):
    results = {}
    kf = KFold(n_splits = settings["kFoldSplits"])

    for train_i, test_i in kf.split(data["x"]):
        testStr = str(test_i[0])+" - " +str(test_i[-1])
        settings["testDataStr"] = testStr

        trainStr = str(train_i[0])+" - "+str(train_i[-1])
        settings["trainDataStr"] = trainStr
        results["powerData"] = data["power_val"]

        x_train, x_test  = data["x"][train_i], data["x"][test_i]
        results["testData"] = x_test
        results["trainData"] = x_train
        if "if" in settings["algo"]:
            outliers = calcIsoForest(x_train, x_test)
            if outliers is None:
                return
            results["outIsoForTest"] = outliers

            outliers = calcIsoForest(x_train, data["power_val"])
            if outliers is None:
                return
            results["outIsoForPower"] = outliers

        if "lof" in settings["algo"]:
            outliers = calcLOF(x_train, x_test)
            if outliers is None:
                return
            results["outLOFTest"] = outliers

            outliers = calcLOF(x_train, data["power_val"])
            if outliers is None:
                return
            results["outLOFPower"] = outliers

        metrics = calcMetrics(results, settings)
        storeResults(settings, metrics)

## GET ALL APPLICATIONS ###

def getApplicationsFromDataset(settings, applicationTypes=None):
    applications = []

    data = DataSet(settings["h5File"])

    for building in data.buildings:
            for appl in data.buildings[building].elec.appliances:
                if (building, appl.type["type"]) not in applications:
                    if applicationTypes is None:
                        applications.append((building,
                                            appl.type["type"]))
                    else:
                        if appl.type["type"] in applicationTypes:
                            applications.append((building,
                                                appl.type["type"]))
    return applications

## GET AGGREGATES ###

def getAggregateFromDataset(settings, applicationTypes=None):
    applications = []

    data = DataSet(settings["h5File"])

    for building in data.buildings:
        if (building, None) not in applications:
            applications.append((building, None))

    return applications

## MAIN ###

if __name__ == "__main__":

    ## sys.argv usage - python outlierMain.py
    #                       h5File
    #                       applications
    #                       windowStart
    #                       windowEnd
    ## applications - ["kettle", "computer"] or
    #                       "aggregate" if aggregate data
    #                       should be used
    ## windowStart - "2013-09-08" or None if whole dataset
    #                       should be used

    if len(sys.argv) != 5:
        print("Not enough arguments")
        sys.exit()

    settings = {}
    files = [sys.argv[1]]

    filters = [None, "rollingMedian", "hampel"]
    applicationTypes = sys.argv[2].split(",")

    for filterType in filters: # None, rollingMedian or hampel
        for h5File in files:
            ### SETTINGS ###
            settings["kFoldSplits"] = 5

            if sys.argv[3] == None:
                settings["window"]      = None
            else:
                settings["window"]      = (sys.argv[3], sys.argv[4])

            settings["filterData"]  = filterType
            settings["algo"]        = "if|lof"
            settings["jsonIndent"]  = False # increases readability
            ### /SETTINGS ###

            ## check if directory with h5Files exist else path passed from user
            if os.path.exists("h5Files"):
                settings["h5File"]      = "h5Files/" + h5File
            else:
                settings["h5File"]      = h5File

            ## create directory for results
            if not os.path.exists("simuResults"):
                os.mkdir("simuResults")

            time = datetime.datetime.now().strftime("_%d%m%Y_%H%M%S.json")
            settings["fileName"]    = "simuResults/simRes_" + h5File + time

            if applicationTypes == "aggregate":
                applications = getAggregateFromDataset(settings)
            else:
                applications = getApplicationsFromDataset(settings, applicationTypes)

            for appl in applications:
                if applicationTypes == "aggregate":
                    data = loadAggregateDataFromDataset(appl, settings)
                else:
                    data = loadApplicationDataFromDataset(appl, settings)

                if data is not None:
                    settings["application"] = appl
                    settings["maxPowVal"]   = int(data["power_val"][-1])
                    calcOutliers(data, settings)


text:
import numpy as np
from numba import njit
from scipy.special import factorial
import timeit

def cumulative_sum_vectorised_math(a, b, c, d, e, f):
    '''
    This function vectorises the mathematical expression in the nested sums
    '''
    exp_min = 5 - (a + c + d + e)
    exp_max = c + d + f + 1
    exp = 2. ** np.arange(exp_min, exp_max)

    b_fac_grid = np.arange(e)
    g_fac_grid = np.arange(a + b)

    fact_b = factorial(b_fac_grid)
    fact_g = factorial(g_fac_grid)

    return cumulative_sum_JIT(a, b, c, d, e, f, exp_min, exp, fact_b, fact_g)

@njit()
def cumulative_sum_JIT(a, b, c, d, e, f, exp_min, exp, fact_b, fact_g):
    Total = 0
    for ai in range(0, a):
        for bi in range(0, b):
            for ci in range(0, c):
                for di in range(0, d):
                    for ei in range(0, e):
                        for fi in range(0, f):
                            for gi in range(0, ai + bi + 1):
                                for hi in range(0, ci + di + 1):
                                    Total += exp[hi - ei + fi - ai - ci - di + 1 - exp_min] * (ei * ei - 2 * (ei * fi) - 7 * di) * fact_b[bi] * fact_g[gi]

    return Total

dim = 5
dim = dim + 1
a = b = c = d = e = f = dim

print(cumulative_sum_vectorised_math(a, b, c, d, e, f))

max_stars_repo_name: adambaskerville/jekyll-theme-chirpy

max_issues_repo_name: adambaskerville/jekyll-theme-chirpy

max_forks_repo_name: adambaskerville/jekyll-theme-chirpy

reason: import numpy,from scipy,from numba


text:
#!/usr/bin/env python3

# This program outputs directory statistics about the current directory.
# Output diagrams of file size distribution and file type distribution are
# generated.
# Steps:
# 1. read current directory contents
# 2. traverse through directory content recursively:
#    if a file is found: store the following info in a list:
#       - full file path with file name
#       - file name
#       - file size
#       - file extension
# 3. aggregate the collected data
# 4. generate figure of:
#   - file name extension distribution histogram of the 20 most frequent file
#     name extensions

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

def process_folder(path):
    current_stats = []
    if os.access(path, os.R_OK) and os.path.isdir(path):
        current_dir = os.listdir(path)
        for file_name in current_dir:
            full_file_path = os.path.join(path, file_name)
            if os.path.isfile(full_file_path):
                file_size = os.stat(full_file_path).st_size
                file_ext = file_name.split(".")[-1]
                current_stat = {
                    "full_file_path": full_file_path,
                    "file_name":      file_name,
                    "file_ext":       file_ext,
                    "file_size":      file_size
                }
                current_stats.append(current_stat)
            else:
                folder_stats = process_folder(full_file_path)
                current_stats += folder_stats
    return current_stats

if __name__ == "__main__":
    path = "./" if len(sys.argv) < 2 else sys.argv[1]
    stats = process_folder(path)
    file_ext_rows_array = [["file_ext", "file_ext_number"]]
    for stat in stats:
        file_ext_rows_array.append([stat["file_ext"], 1.0])

    file_ext_rows_np_array = np.array(file_ext_rows_array)

    file_ext_pd_dataframe = pd.DataFrame(
        data    = file_ext_rows_np_array[1:,0:],
        columns = file_ext_rows_np_array[0,0:])

    file_ext_pd_dataframe.sort_values(
        by = ["file_ext"],
        ascending = True)

    file_ext_pd_aggregation = {
        "file_ext": {
            "number_of_files": "sum"
        }
    }

    file_ext_count_pd_dataframe = file_ext_pd_dataframe.groupby("file_ext").count()

    file_ext_count_pd_dataframe = file_ext_count_pd_dataframe.sort_values(
        by = ["file_ext_number"],
        ascending = False)

    file_ext_head = file_ext_count_pd_dataframe.head(20)
    file_ext_numbers_head = file_ext_head['file_ext_number'].tolist()
    file_ext_names_head = file_ext_head.index.tolist()

    barplot = plt.bar(range(len(file_ext_head)), file_ext_numbers_head, color='#4CC4FF')
    plt.ylim([0, max(file_ext_numbers_head)])
    xticks = plt.xticks(range(len(file_ext_head)), file_ext_names_head, rotation=45, color="#2BD5E1")
    plt.savefig('dirstats.png', bbox_inches='tight')

max_stars_repo_name: nlharri/dirstats-python

max_issues_repo_name: nlharri/dirstats-python

max_forks_repo_name: nlharri/dirstats-python


text:
#! /usr/bin/python3

import os
import time
import numpy as np
import math as m
import matplotlib.pyplot as plt
from   matplotlib import animation

from xarray import open_mfdataset

# mesonhFiles = '/home/pnarvor/work/nephelae/data/MesoNH-2019-02/REFHR.1.ARMCu.4D.nc'
# mesonhFiles = '/home/pnarvor/work/nephelae/data/nephelae-remote/MesoNH02/bomex_hf.nc'
mesonhFiles = '/home/pnarvor/work/nephelae/data/skyscanner-remote/Nephelae/MesoNH-2019-10/L12zo.1.BOMEX.OUT.040.nc'

atm = open_mfdataset(mesonhFiles)
wt  = atm.WT
print("Shape :", wt.shape)

# t0 = 40
t0 = 0
z0 = 40

borderSize = 10
snBorders  = np.empty([wt.shape[2], 2*borderSize])
snBorders[:, :borderSize] = wt[t0,z0,:,-borderSize:]
snBorders[:,-borderSize:] = wt[t0,z0,:, :borderSize]

toRemove = 6
wt0 = wt[t0,z0,:,:-toRemove]
snBorders0  = np.empty([wt0.shape[0], 2*borderSize])
snBorders0[:, :borderSize] = wt0[:,-borderSize:]
snBorders0[:,-borderSize:] = wt0[:, :borderSize]

pad = 0
bottom = np.array(wt[t0,z0,:,:(toRemove + pad)]).squeeze()
top    = np.array(wt[t0,z0,:,-(toRemove + pad):]).squeeze()
diff0   = top - bottom
print("Residual top/bottom:", np.linalg.norm(diff0.ravel()))

pad = 0
left  = np.array(wt[t0,z0,:(toRemove + pad) ,:]).squeeze()
right = np.array(wt[t0,z0,-(toRemove + pad):,:]).squeeze()
diff1 = right - left
print("Residual right/left:", np.linalg.norm(diff1.ravel()))

fig, axes = plt.subplots(3,1, sharex=True)
axes[0].imshow(np.array(wt[t0,z0,:,:].data).squeeze().T, origin='lower')
axes[1].imshow( snBorders.T, origin='lower')
axes[2].imshow(snBorders0.T, origin='lower')

fig, axes = plt.subplots(3,1, sharex=True)
axes[0].set_title("Bottom/Top")
axes[0].imshow(bottom.T, origin='lower')
axes[1].imshow(   top.T, origin='lower')
axes[2].imshow( diff0.T, origin='lower')

fig, axes = plt.subplots(3,1, sharex=True)
axes[0].set_title("Left/Right")
axes[0].imshow( left, origin='lower')
axes[1].imshow(right, origin='lower')
axes[2].imshow(diff1, origin='lower')

plt.show(block=False)


text:
from __future__ import division

import math

import numpy as np

# Non-monotonic Ishigami Function (3 parameters)
# First-order indices:
# x1: 0.3139
# x2: 0.4424
# x3: 0.0
def evaluate(values):
    Y = np.zeros([values.shape[0]])
    A = 7
    B = 0.1

    for i, X in enumerate(values):
        Y[i] = math.sin(X[0]) + A * math.pow(math.sin(X[1]), 2) + \
            B * math.pow(X[2], 4) * math.sin(X[0])

    return Y


text:
#!/usr/bin/env python

"""
How many numbers are there between 1 and 1,000,000 whose sum of its individual
digits are equal to the last two digits?

Example: 512119 qualifies because 5 + 1 + 2 + 1 + 1 + 9 == 19

This is a literal translation of "Number Problem.cpp" in this same repo, which
I wrote a long time ago for some reason.  I've changed the variable naming
convention from CamelCase to underscores for sanity.
"""

from numba import autojit, jit

@jit("void()")
#@autojit()
def compute():
    numbercount = 0

    for i in range(1, 1000001):
        ones = i % 10
        tens = i % 100 // 10
        hundreds = i % 1000 // 100
        thousands = i % 10000 // 1000
        ten_thousands = i % 100000 // 10000
        hundred_thousands = i % 1000000 // 100000
        millions = i % 10000000 // 1000000

        digit_sum = (ones + tens + hundreds + thousands + ten_thousands + hundreds +
                     hundred_thousands + millions)

        last_two_digits = tens*10 + ones

        if digit_sum == last_two_digits:
            #print("%d works with %d!" % (i, digit_sum))
            numbercount += 1

    print("The number of numbers is %d." % (numbercount,))

for i in range(10):
    compute()

max_stars_repo_path: Number Problem-numba.py

max_stars_repo_name: asmeurer/Miscellaneous-Useful-Scripts

max_issues_repo_path: Number Problem-numba.py

max_issues_repo_name: asmeurer/Miscellaneous-Useful-Scripts

max_forks_repo_path: Number Problem-numba.py

max_forks_repo_name: asmeurer/Miscellaneous-Useful-Scripts

reason: from numba


text:
import tensorflow as tf
import numpy as np
import datasets
from utils import config, Struct
import utils
import os
import sys
from nnet import Affine, Chain, LeakyRelu, BernoulliLogits, Const, ZeroBias
import matplotlib
matplotlib.use('Agg')
import matplotlib.pyplot as plt

class Model():
    def __init__(self, x, latent_units, q_units, p_units, sleep='usual', batch_norm=True):
        self.x = x
        self.mb_size = tf.shape(x)[0]
        self.latent_units = latent_units
        datadim = self.x.get_shape().as_list()[1]
        self.q_units = []
        for n1, n2, u in zip([datadim] + latent_units, latent_units, q_units):
            self.q_units.append([n1] + u + [n2])
        self.p_units = []
        for n1, n2, u in zip(reversed(latent_units), list(reversed(latent_units))[1:] + [datadim], p_units):
            self.p_units.append([n1] + u + [n2])
        self._setup_connectivity(batch_norm)
        self._wake_phase()
        self._sleep_phase(sleep)

        self.loss = - (self.log_q_sleep + self.log_p_wake)

        self._test_loss()

    def _setup_connectivity(self, batch_norm):
        self.alpha = tf.Variable(0.)
        self.q = []
        for units in self.q_units:
            self.q.append(
                Chain(Chain.build_chain(units, LeakyRelu, batch_norm=batch_norm, batch_norm_alpha=self.alpha),
                      BernoulliLogits()
                      )
            )
        self.p = []
        for units in self.p_units:
            self.p.append(
                Chain(Chain.build_chain(units, LeakyRelu, batch_norm=batch_norm, batch_norm_alpha=self.alpha),
                      BernoulliLogits()
                      )
            )
        self.prior = Chain(Const(ZeroBias(self.latent_units[-1])),
                           BernoulliLogits()
                           )

    def _wake_phase(self):
        self.q_samples = [self.x]
        for q in self.q:
            self.q_samples.append(tf.stop_gradient(q(self.q_samples[-1]).out.sample))
        self.log_p_wake = self.prior(self.mb_size).out.log_prob(self.q_samples[-1])
        for p, sample1, sample2 in zip(reversed(self.p), self.q_samples[1:], self.q_samples):
            self.log_p_wake += p(sample1).out.log_prob(sample2)

    def _sleep_phase(self, sleep_kind):
        # if not (sleep_kind == 'top_layer_mix' or sleep_kind == 'top_down_mix'):
        self.prior_samples = [self.prior(self.mb_size).out.sample]
        for p in self.p:
            self.prior_samples.append(p(self.prior_samples[-1]).out.sample)

        # this is the crucial part
        if sleep_kind == 'usual':
            self.sleep_samples = self.prior_samples
        elif sleep_kind == 'autoenc':
            self.sleep_samples = list(reversed(self.q_samples))
        elif sleep_kind == 'mix':
            # when alpha is zero, take the sample from q; when 1 - from the prior
            self.mask = tf.cast(tf.random_uniform([1]) < self.alpha, tf.float32)
            self.sleep_samples = [self.mask * prior_sample + (1 - self.mask) * q_sample
                                  for prior_sample, q_sample in zip(self.prior_samples, reversed(self.q_samples))]
        elif sleep_kind == 'mixy_mix':
            # when alpha is zero, take the sample from q; when 1 - from the prior
            self.mask = tf.cast(tf.random_uniform([len(self.latent_units)+1,1]) < self.alpha, tf.float32)
            self.sleep_samples = [mask * prior_sample + (1 - mask) * q_sample
                                  for mask, prior_sample, q_sample in zip(self.mask, self.prior_samples, reversed(self.q_samples))]
        elif sleep_kind == 'reconstruct':
            # p samples sleep_samples from last q sample instead of prior sample
            self.sleep_samples = [self.q_samples[-1]]
            for p in self.p:
                self.sleep_samples.append(p(self.sleep_samples[-1]).out.sample)
        elif sleep_kind == 'top_down_mix1':
            # sample sleep from q sample at floor(alpha*10) layer
            layer = len(self.latent_units) - tf.floor((len(self.latent_units)+1)*self.alpha)
            i = tf.constant(0.)
            self.mask = tf.cast((layer < i),tf.float32)
            self.sleep_samples = [(1-self.mask) * self.q_samples[-1] + self.mask * self.prior_samples[0]]
            for q_sample, p in zip(reversed(self.q_samples[:-1]), self.p):
                self.sleep_samples.append( self.mask * p(self.sleep_samples[-1]).out.sample + (1-self.mask) * q_sample)
                i = i+1
                self.mask = tf.cast((layer < i),tf.float32)
        elif sleep_kind == 'top_down_mix2':
            # use alpha to determine prior sample or q samples. if q samples then sample sleep from q sample at floor(alpha*10) layer
            if tf.cast(tf.random_uniform([1]) < self.alpha, tf.int32)==1.: self.sleep_samples = self.prior_samples
            else:
                layer = len(self.latent_units) - tf.floor((len(self.latent_units)+1)*self.alpha)
                i = tf.constant(0.)
                self.mask = tf.cast((layer < i),tf.float32)
                self.sleep_samples = [(1-self.mask) * self.q_samples[-1] + self.mask * self.prior_samples[0]]
                for q_sample, p in zip(reversed(self.q_samples[:-1]), self.p):
                    self.sleep_samples.append( self.mask * p(self.sleep_samples[-1]).out.sample + (1-self.mask) * q_sample)
                    i = i+1
                    self.mask = tf.cast((layer < i),tf.float32)

        self.log_q_sleep = 0
        for q, p, sample in zip(reversed(self.q), self.p, self.sleep_samples):
            sample1 = tf.stop_gradient(p(sample).out.sample)
            self.log_q_sleep += q(sample1).out.log_prob(sample)

    def _test_loss(self):
        self.k = tf.placeholder(tf.int32)
        self.x_rep = utils.tf_repeat(self.x, self.k)
        self.q_samples_rep = [self.x_rep]
        for q in self.q:
            self.q_samples_rep.append(q(self.q_samples_rep[-1]).out.sample)

        self.variational_lower_bound = self.prior(self.k*self.mb_size).out.log_prob(self.q_samples_rep[-1])
        for p, sample1, sample2 in zip(reversed(self.p), self.q_samples_rep[1:], self.q_samples_rep):
            self.variational_lower_bound += p(sample1).out.log_prob(sample2)
        for q, sample1, sample2 in zip(self.q, self.q_samples_rep, self.q_samples_rep[1:]):
            self.variational_lower_bound -= q(sample1).out.log_prob(sample2)

        self.variational_lower_bound = tf.reshape(self.variational_lower_bound, tf.pack([-1, self.k]))
        self.variational_lower_bound = utils.tf_log_mean_exp(self.variational_lower_bound)

        self.reconstruction_loss = self.p[-1](self.q_samples_rep[1]).out.log_prob(self.x_rep)
        self.reconstruction_loss = tf.reshape(self.reconstruction_loss, tf.pack([-1, self.k]))
        self.reconstruction_loss = tf.reduce_mean(self.reconstruction_loss, 1)

def main(args):
    args.update_exclusive(default_args)
    directory = os.path.join(config.RESULTSDIR, 'ws', args.directory)
    if not os.path.exists(directory):
        os.makedirs(directory)
    if not args.restore:
        map(os.remove, (os.path.join(directory, f) for f in os.listdir(directory) if f.endswith(".txt") or f.endswith(".png")))
        with open(os.path.join(directory, 'params.txt'), 'w') as f:
            f.write(repr(args))
    print directory

    dataset = datasets.MNIST(binary=True)

    x = tf.placeholder(np.float32, shape=(None, dataset.get_data_dim()))

    model = Model(x, args.latent_units,
                  q_units=args.q_units,
                  p_units=args.p_units,
                  sleep=args.sleep_type,
                  batch_norm=args.bn)
    examples_per_epoch = dataset.data['train'][0].shape[0]
    num_updates = args.n_epochs * examples_per_epoch / args.mb_size

    step = tf.Variable(0, trainable=False)
    lr = tf.placeholder(tf.float32)
    train_op = args.optimizer(lr).minimize(model.loss, global_step=step)
    with tf.control_dependencies([train_op]):
        # linearly anneal alpha from 0 to 1 over the course of N/2 epochs, then train for additional N/2 epochs with alpha=1
        train_op = tf.assign(model.alpha, tf.minimum(1., tf.maximum(0., tf.cast(step, tf.float32) / num_updates * 2)))

    init_op = tf.initialize_all_variables()

    saver = tf.train.Saver(max_to_keep=1)

    with tf.Session() as sess:
        if not args.restore:
            sess.run(init_op)
        else:
            saver.restore(sess, tf.train.latest_checkpoint(os.path.join(directory)))
            print "restored"
            measure_test_log_likelihood(sess, model, dataset, directory)
            import ipdb; ipdb.set_trace()
        for x_np, _ in dataset.random_minibatches('train', args.mb_size, num_updates):
            i, _ = sess.run([step, train_op], feed_dict={x: x_np, lr: args.lr})
            if i % 1000 == 1 or i == num_updates - 1:
                visualize(sess, model, dataset, directory, float(i) * args.mb_size / examples_per_epoch)
            if i % 10000 == 1 or i == num_updates - 1:
                saver.save(sess, os.path.join(directory, 'model.chk'), global_step=step)
            if i % 10000 == 1:
                print directory

def visualize(sess, model, dataset, directory, epoch):
    x_np, _ = dataset.get_random_minibatch('train', 100, np.random.RandomState(123))
    x_np_valid, _ = dataset.get_random_minibatch('valid', 100, np.random.RandomState(123))
    q_means = [q(s).out.mu for q, s in zip(model.q, model.q_samples)]
    q_means = sess.run(q_means, feed_dict={model.x: x_np})
    for i, q_m in enumerate(q_means):
        plt.boxplot(q_m[:, np.argsort(q_m.std(axis=0))])
        # plt.violinplot(q_means)
        plt.savefig(os.path.join(directory, 'means of q{}.png'.format(i)))
        plt.close()
    prior_means = sess.run(model.prior(1).out.mu)[0]
    plt.scatter(np.arange(prior_means.shape[0]), prior_means)
    plt.savefig(os.path.join(directory, 'prior means.png'))
    plt.close()
    # pre_activations = [layer.y for layer in model.q_params_chain.all_layers]
    # names = [layer.to_short_str() for layer in model.q_params_chain.all_layers]
    # pre_activations = sess.run(pre_activations, feed_dict={model.x: x_np})
    # for i, (act, name) in enumerate(zip(pre_activations, names)):
    #     plt.violinplot(act[:, np.argsort(act.std(axis=0))])
    #     plt.savefig(os.path.join(directory, 'rec_activations{}{}.png'.format(i+1, name)))
    #     plt.close()
    # pre_activations = [layer.y for layer in model.p_samples_copy.layers]
    # names = [layer.to_short_str() for layer in model.p_samples_copy.layers]
    # pre_activations = sess.run(pre_activations, feed_dict={model.num_samples: 100})
    # for i, (act, name) in enumerate(zip(pre_activations, names)):
    #     plt.violinplot(act[:, np.argsort(act.mean(axis=0))])
    #     plt.savefig(os.path.join(directory, 'gen_model_activations{}{}.png'.format(i+1, name)))
    #     plt.close()
    first_linear_layer = next(layer for layer in model.q[0].layers[0].layers if isinstance(layer, Affine))
    order = np.argsort(sess.run(model.q[0].layers[0].layers[0](model.x).y, feed_dict={model.x: x_np}).std(axis=0))
    W = sess.run(first_linear_layer.W)
    plt.imshow(utils.misc.tile_images(dataset.reshape_for_display(W.T[order])), cmap='gray')
    plt.savefig(os.path.join(directory, 'rec_weights.png'))
    plt.close()
    last_linear_layer = [layer for layer in model.p[-1].layers[0].layers if isinstance(layer, Affine)][-1]
    W = sess.run(last_linear_layer.W)
    plt.imshow(utils.misc.tile_images(dataset.reshape_for_display(W)), cmap='gray')
    plt.savefig(os.path.join(directory, 'gen_weights.png'))
    plt.close()
    samples, recs = sess.run([model.p[-1](model.prior_samples[-2]).out.mu,
                              model.p[-1](model.q[0](model.x).out.sample).out.mu], feed_dict={model.x: x_np})
    samples = utils.misc.tile_images(dataset.reshape_for_display(samples))
    recs = utils.misc.tile_images(dataset.reshape_for_display(recs))
    orig = utils.misc.tile_images(dataset.reshape_for_display(x_np))
    orig_and_samples = np.concatenate([orig, recs, samples], axis=1)
    plt.imshow(orig_and_samples, cmap='gray')
    plt.savefig(os.path.join(directory, 'samples.png'))
    plt.close()

    loss_1 = sess.run(model.variational_lower_bound, feed_dict={model.x: x_np, model.k: 1})
    loss_10 = sess.run(model.variational_lower_bound, feed_dict={model.x: x_np, model.k: 10})
    loss_100, rec_100 = sess.run([model.variational_lower_bound, model.reconstruction_loss],
                                 feed_dict={model.x: x_np, model.k: 100})
    loss_valid, rec_valid = sess.run([model.variational_lower_bound, model.reconstruction_loss],
                                     feed_dict={model.x: x_np_valid, model.k: 100})
    utils.print_and_save(os.path.join(directory, 'loss.txt'),
                         'epoch ', epoch,
                         'train rec100', np.mean(rec_100),
                         'valid rec100', np.mean(rec_valid),
                         'train loss1', np.mean(loss_1),
                         'train loss10', np.mean(loss_10),
                         'train loss100', np.mean(loss_100),
                         'valid loss100', np.mean(loss_valid),
                         'alpha', sess.run(model.alpha))

    # loss = measure_validation_loss(sess, model, dataset)
    # utils.print_and_save(os.path.join(directory, 'loss.txt'), 'train', cum_loss, 'valid', loss)

def measure_test_log_likelihood(sess, model, dataset, directory, k=5000):
    mbsize = 10
    ll = 0.
    n = 0
    for x_np, _ in dataset.all_minibatches('test', mbsize, np.random.RandomState(123)):
        n += x_np.shape[0]
        ll += np.sum(sess.run(model.variational_lower_bound, feed_dict={model.x: x_np, model.k: k}))
        sys.stdout.write(str(n)+','+str(ll/n))
    utils.print_and_save(os.path.join(directory, 'loss.txt'),
                         'test log likelihood 5000', ll/n)

default_args = Struct(
    optimizer=tf.train.AdamOptimizer,
    sleep_type='usual',
    mb_size=200,
    latent_units=50,
    q_units=[500, 500],
    p_units=[500, 500],
    bn=True,
    lr=0.001,
    n_epochs=500,
    directory='default',
    restore=False
    )

max_stars_repo_name: jackklys/wake-sleep-modifications

max_issues_repo_name: jackklys/wake-sleep-modifications

max_forks_repo_name: jackklys/wake-sleep-modifications


text:
import geopandas as gpd
import numpy as np
from libpysal.examples import load_example
from segregation.batch import (
    batch_compute_multigroup,
    batch_compute_singlegroup,
    batch_multiscalar_singlegroup,
    batch_multiscalar_multigroup,
)

s_map = gpd.read_file(load_example("Sacramento1").get_path("sacramentot2.shp"))

def test_batch_single():
    fit = batch_compute_singlegroup(
        s_map.to_crs(s_map.estimate_utm_crs()),
        group_pop_var="HISP",
        total_pop_var="TOT_POP",
        distance=2000,
        center="mean",
        function="triangular",
    )
    np.testing.assert_array_almost_equal(
        fit.Statistic,
        [
            0.6905,
            0.061,
            0.8511,
            0.1218,
            0.2922,
            0.2639,
            0.1243,
            0.0746,
            0.8043,
            0.22,
            0.2924,
            0.86,
            0.2,
            0.0767,
            0.3932,
            0.7717,
            0.2283,
            0.4524,
            0.2839,
            0.3827,
            0.3113,
            -0.1034,
            0.52,
            0.1264,
            0.2612,
            0.2285,
            1.0468,
        ],
        decimal=2,
    )

def test_batch_multi():
    mfit = batch_compute_multigroup(
        s_map.to_crs(s_map.estimate_utm_crs()),
        distance=2000,
        groups=["HISP", "BLACK", "WHITE"],
    )
    np.testing.assert_array_almost_equal(
        mfit.Statistic,
        [
            0.37768411,
            0.11294892,
            0.78242435,
            0.50485431,
            0.14435762,
            0.15922567,
            0.13971496,
            0.11631982,
            0.55281779,
            0.44718221,
        ],
        decimal=3,
    )

def test_batch_multiscalar_multi():
    mfit = batch_multiscalar_multigroup(
        s_map.to_crs(s_map.estimate_utm_crs()),
        distances=[500, 1000],
        groups=["HISP", "BLACK", "WHITE"],
    )
    assert mfit.shape == (3, 10)

def test_batch_multiscalar_single():
    mfit = batch_multiscalar_singlegroup(
        s_map.to_crs(s_map.estimate_utm_crs()),
        distances=[500, 1000],
        group_pop_var="HISP",
        total_pop_var="TOT_POP",
    )
    assert mfit.shape == (3, 13)

max_stars_repo_name: noahbouchier/segregation

max_issues_repo_name: noahbouchier/segregation

max_forks_repo_name: noahbouchier/segregation


text:
import matplotlib.colors
import matplotlib.pyplot as plt
import numpy as np
from pcl_helper import *

def rgb_to_hsv(rgb_list):
    rgb_normalized = [1.0*rgb_list[0]/255, 1.0*rgb_list[1]/255, 1.0*rgb_list[2]/255]
    hsv_normalized = matplotlib.colors.rgb_to_hsv([[rgb_normalized]])[0][0]
    return hsv_normalized

### utility to calculate normailized histogram
##  params: 3 channel values
##          nbins: bins size in histogram, default to 32
##	    bins_range: default to (0, 256)
def get_normalized_hist(channel_1_vals, channel_2_vals, channel_3_vals, nbins=32, bins_range=(0,256)):
    # Compute histograms
    hist0 = np.histogram(channel_1_vals, bins=nbins, range=bins_range)
    hist1 = np.histogram(channel_2_vals, bins=nbins, range=bins_range)
    hist2 = np.histogram(channel_3_vals, bins=nbins, range=bins_range)

    # Concatenate and normalize the histograms
    hist_features = np.concatenate((hist0[0], hist1[0], hist2[0])).astype(np.float64)

    return  hist_features / np.sum(hist_features)

def compute_color_histograms(cloud, nbins=32, using_hsv=False):

    # Compute histograms for the clusters
    point_colors_list = []

    # Step through each point in the point cloud
    for point in pc2.read_points(cloud, skip_nans=True):
        rgb_list = float_to_rgb(point[3])
        if using_hsv:
            point_colors_list.append(rgb_to_hsv(rgb_list) * 255)
        else:
            point_colors_list.append(rgb_list)

    # Populate lists with color values
    channel_1_vals = []
    channel_2_vals = []
    channel_3_vals = []

    for color in point_colors_list:
        channel_1_vals.append(color[0])
        channel_2_vals.append(color[1])
        channel_3_vals.append(color[2])

    normed_features = get_normalized_hist(channel_1_vals, channel_2_vals, channel_3_vals, nbins=nbins)

    return normed_features

def compute_normal_histograms(normal_cloud, nbins=32):

    # Populate lists with normal values
    norm_x_vals = []
    norm_y_vals = []
    norm_z_vals = []

    for norm_component in pc2.read_points(normal_cloud,
                                          field_names = ('normal_x', 'normal_y', 'normal_z'),
                                          skip_nans=True):
        norm_x_vals.append(norm_component[0])
        norm_y_vals.append(norm_component[1])
        norm_z_vals.append(norm_component[2])

    # get normalized features with bins_range = (-1,1)
    normed_features = get_normalized_hist(norm_x_vals, norm_y_vals, norm_z_vals, nbins=nbins, bins_range=(-1,1))

    return normed_features


text:
# -*- coding: utf-8 -*-
# Based on https://stackoverflow.com/questions/24633664/confidence-interval-for-exponential-curve-fit/37080916#37080916

import matplotlib.pyplot as plt
import numpy as np
from datetime import date
from scipy.stats import describe
import os

def figurize(mdl,sfx,guess=None):
    # Arguments are a model object and a suffix for the filename.
    fig = plt.figure(figsize=(6, 4))
    plt.suptitle("COVID-19 Cases: District of Columbia", fontweight="bold")
    plt.title("github.com/reidac/covid19-curve-dc", style="oblique")
    plt.xlabel("Days since March 8, 2020")
    plt.ylabel("# Diagnosed Cases")

    # Data

    import get_dc_data
    casedata = get_dc_data.retrieve()

    # The new improved CaseData object promises us data
    # sorted in order of increasing x.
    x = casedata.x
    y = casedata.positive

    plt.scatter(x, y, marker=".", s=10, color="k", zorder=10)

    # Levenburg-Marquardt Least-Squares Fit

    if (guess is not None):
        mdl.params = guess
    mdl.fit(x,y)

    # print("Params: ", mdl.params)

    # print("cases ~ {0:.2g} * (1 + {1:.2g})^t".format(a, b))

    # Confidence Band: dfdp represents the partial derivatives of
    # the model with respect to each parameter p (i.e., a and b)

    # Use casedata, not the trailing item in the list, which for us
    # is not the largest abscissa!
    xhat = np.linspace(0, casedata.today-casedata.start+7, 100)
    yhat = mdl.evaluate(xhat)

    upper = mdl.evaluate_h(xhat)
    lower = mdl.evaluate_l(xhat)

    ix = np.argsort(xhat)
    plt.plot(xhat[ix], yhat[ix], c="#FFCCCC", lw=1, zorder=5)
    plt.fill_between(
        xhat[ix], upper[ix], yhat[ix], edgecolor=None,
        facecolor="#F0F0F0", zorder=1
    )
    plt.fill_between(
        xhat[ix], lower[ix], yhat[ix], edgecolor=None,
        facecolor="#F0F0F0", zorder=1
    )

    # Plot Boundaries

    plt.xlim([0, xhat[-1]])
    plt.ylim([0, max(upper[-1]*1.1,y[-1]*1.1)])

    xmax = xhat[-1]
    ymax = upper[-1]

    # Predictions

    tomorrow = date.fromordinal(casedata.today + 1)
    nextWeek = date.fromordinal(casedata.today + 7)

    # print(tomorrow.toordinal()-casedata.start)
    # print(nextWeek.toordinal()-casedata.start)
    xhat = np.array([tomorrow.toordinal() - casedata.start,
                     nextWeek.toordinal() - casedata.start])

    yhat = mdl.evaluate(xhat)
    upper = mdl.evaluate_h(xhat)
    lower = mdl.evaluate_l(xhat)

    dx = 0.25

    plt.text(
        dx,
        2.0*ymax/3.0,
        "{0}/{1}: ({2:.0f} < {3:.0f} < {4:.0f})".format(
            tomorrow.month, tomorrow.day, lower[0], yhat[0], upper[0]
        ),
        va="center",
        zorder=5,
        bbox=dict(boxstyle="round", ec="black", fc="white", linewidth=dx),
    )
    plt.text(
        dx,
        2.5*ymax/3.0,
        "{0}/{1}: ({2:.0f} < {3:.0f} < {4:.0f})".format(
            nextWeek.month, nextWeek.day, lower[1], yhat[1], upper[1]
        ),
        va="center",
        zorder=5,
        bbox=dict(boxstyle="round", ec="black", fc="white", linewidth=dx),
    )

    plt.arrow(
        0.4*xmax,
        (2.0*ymax/3.0),
        xhat[0] - 0.0625 - 0.4*xmax,
        yhat[0]-(2.0*ymax/3.0),
        fc="black",
        ec="black",
        length_includes_head=True,
        linewidth=0.5,
        linestyle="--",
        zorder=2,
    )
    plt.arrow(
        0.4*xmax,
        2.5*(ymax/3.0),
        xhat[1] - 0.0625 - 0.4*xmax,
        yhat[1]-(2.5*ymax/3.0),
        fc="black",
        ec="black",
        length_includes_head=True,
        linewidth=0.5,
        linestyle="--",
        zorder=2,
    )

    # Save figure

    if "FIG_PATH" in os.environ:
        fig_path = os.environ['FIG_PATH']
    else:
        fig_path = "."

        plt.savefig("{0}/us_dc_{1}.png".format(fig_path,sfx), dpi=300,
                    bbox_inches="tight")

if __name__=="__main__":

    import model
    m = model.Exponential()
    figurize(m,"exp")

reason: import numpy,from scipy


text:
# -*- coding: utf-8 -*-
"""
Created on Thu Nov 18 15:43:55 2021

@author: ZeitgeberH
"""
from pathlib import Path
from PyQt5 import QtGui, QtCore, QtSvg
from pyqtgraph.Qt import QtWidgets
from PyQt5.QtWidgets import QMessageBox, QTableWidgetItem
import pyqtgraph as pg
import pyqtgraph.opengl as gl
from pyqtgraph.parametertree import Parameter, ParameterTree
import gr_pars as AllMyPars
from pyqtgraph import GraphicsLayoutWidget
import sys
import os
import numpy as np
import pandas as pd
import itertools
import time as sysTime
from copy import deepcopy
import glob
from gr_viewClass import (FileView,ParameterView,PlotView,
    TabView,TreeView, TableView,showdialog,getfiles,ImageView_gv)
from gr_allenFISH import GeneRenderAPI
from gr_allenFISH_utils import (
    check_gene_cached,
    check_gene_cached_images,
    load_cached_gene,
    download_and_cache,
    download_and_cache_image,
    queryGene_ncbi,rgb2intensity
)
from gr_sys_utils import base_dir
import pdb
from PIL import Image
Image.MAX_IMAGE_PIXELS = 1000000000  ## silence decompression bomb warning
import tqdm
import webbrowser
from skimage import filters as SKIfilters
from skimage import measure as SKImeasure
from matplotlib.cm import jet as jetcm

class MainWindow(QtWidgets.QMainWindow):
    """
    Main frame.
    """
    def __init__(self, app,  parent = None):
        super(MainWindow, self).__init__(parent)
        self.app = app
        self.create_mainWindow()
        self.setWindowTitle("FISHViewer")
        self.setWindowIcon(pg.QtGui.QIcon('data\\icons\\Fish.png'))
        self.setWindowState(QtCore.Qt.WindowMaximized)

    def create_mainWindow(self):
        # Configure Qt GUI:
        self.make_layout()
        self.add_menubar()
        self.add_toolbar()
        self.statusBar = QtWidgets.QStatusBar()
        self.setStatusBar(self.statusBar)
        self.MouseMode =  pg.ViewBox.RectMode
        self.mainFrame = pg.QtGui.QWidget()
        self.mainFrame.setLayout(self.MainLayout)
        self.setCentralWidget(self.mainFrame)
        self.otherInitStuff()
        self.setRightPanelStrechFactor()

    def make_layout(self):
        ## layout
        self.MainLayout = pg.QtGui.QGridLayout()
        # horizontal splitter
        self.frame_splitter = QtWidgets.QSplitter(QtCore.Qt.Horizontal)
        # splitter for file browswer (top) and  pul view (bottom)
        self.right_panels = QtWidgets.QSplitter(QtCore.Qt.Vertical, parent = self.frame_splitter)

        ### left pannels for main plotting and data view
        self.visulization_view = TabView(self.frame_splitter)
        self.plot_splitter = QtWidgets.QSplitter(QtCore.Qt.Vertical, parent = self.visulization_view)

        self.ISH_view = pg.ImageView(view=pg.PlotItem()) #, SynAtlasFunc = self.SynAtlas)
        self.ISH_view_vLine = pg.InfiniteLine(angle=90, movable=False, name='cursorV1')
        self.ISH_view_hLine = pg.InfiniteLine(angle=0, movable=False, name='cursorH1')
        self.ISH_view.addItem(self.ISH_view_vLine, ignoreBounds=True)
        self.ISH_view.addItem(self.ISH_view_hLine, ignoreBounds=True)
        self.ISH_view.view.hoverEvent = self.imageHoverEvent

        self.ISH_view.view.scene().sigMouseClicked.connect(self.ISH_view_mouseClicked)
        # self.ISH_view.ui.histogram.hide()

        self.plot_splitter.addWidget(self.ISH_view)
        # self.updateFishImage(self.ISH_view, 'atlas_example.jpg')
        self.ISH_view.setColorMap(pg.colormap.get('CET-L1'))
        self.ISH_view.ui.roiBtn.hide()
        self.ISH_view.ui.menuBtn.hide()
        self.ISH_view.roi.hide()
        # self.ISH_view.ui.roiPlot.close()
        # self.ISH_view.ui.roiPlot.hide()
        self.scale1 = pg.ScaleBar(size=50, suffix='pixels')
        self.scale1.setParentItem(self.ISH_view.view.getViewBox())
        self.scale1.anchor((1, 1), (1, 1), offset=(-20, -20))

        self.Expression_view = pg.ImageView(view=pg.PlotItem())
        self.plot_splitter.addWidget(self.Expression_view)
        self.Expression_view_vLine = pg.InfiniteLine(angle=90, movable=False, name='cursorV2')
        self.Expression_view_hLine = pg.InfiniteLine(angle=0, movable=False, name='cursorH2')
        self.Expression_view.addItem(self.Expression_view_vLine, ignoreBounds=True)
        self.Expression_view.addItem(self.Expression_view_hLine, ignoreBounds=True)
        self.Expression_view.view.hoverEvent = self.imageHoverEvent
        # self.updateFishImage(self.Expression_view, 'atlas_example.jpg')
        self.Expression_view.setColorMap(pg.colormap.get("CET_L16", source='colorcet')) #pg.colormap.get('CET-L1'))
        self.Expression_view.ui.roiBtn.hide()
        self.Expression_view.ui.menuBtn.hide()
        self.Expression_view.roi.hide()
        self.Expression_view.ui.roiPlot.getPlotItem().hideAxis('bottom')
        self.scale2 = pg.ScaleBar(size=50, suffix='pixels')
        self.scale2.setParentItem(self.Expression_view.view.getViewBox())
        self.scale2.anchor((1, 1), (1, 1), offset=(-20, -20))
        # self.Expression_view.ui.roiPlot.close()
        # self.Expression_view.ui.roiPlot.hide()

        self.Expression_view.view.setXLink(self.ISH_view.view)
        self.Expression_view.view.setYLink(self.ISH_view.view)

        self.visulization_view.addTab(self.plot_splitter, 'ABI Views')

        self.frame_splitter.addWidget(self.visulization_view)

        ## right panels for parameters and input
        self.frame_splitter.addWidget(self.right_panels)
        self.files_view = TabView(self.right_panels)

        ## parameters for ABI FISH images
        self.globalPars = Parameter.create(name='params', type='group', children=AllMyPars.download_options, readonly=False)
        self.globalPars.sigTreeStateChanged.connect(self.event_parameters_stateChange) ## call back
        self.parTree = ParameterTree()
        self.parTree.setHeaderLabels(["Parameter                                  ", "Value"])
        self.parTree.setParameters(self.globalPars, showTop=False)
        self.parTree.setWindowTitle('Global parameter')
        self.files_view.addTab(self.parTree, 'ABA parameters')

        self.right_panels.addWidget(self.files_view)

        ### Middle panel: recording paramters
        self.experiments_view = TabView(self.right_panels)
        self.right_panels.addWidget(self.experiments_view)

        self.experimentsInfo_Tab = TableView(self.experiments_view,editable=False, sortable=False, mdcFunc=self.loadSectionImageInfo)
        self.experiments_view.addTab(self.experimentsInfo_Tab, 'Experiments')

        self.experimentImage_Tab = TableView(self.experiments_view,editable=False, sortable=False) #, mdcFunc=self.render_sectionImage,kpFunc=self.render_sectionImage)
        self.experimentImage_Tab.currentCellChanged.connect(self.render_sectionImage)
        self.experiments_view.addTab(self.experimentImage_Tab, 'Section images')

        ## bottom panel: extra informations
        self.trees_view = TabView(self.right_panels)
        self.geneInfo_Tab = TableView(self.trees_view,mdcFunc=self.querySelectedGeneInfo) ## Table to show queried gene information
        self.trees_view.addTab(self.geneInfo_Tab,'Gene information' )

        self.abAtlas_view = TableView(self.trees_view)
        self.trees_view.addTab(self.abAtlas_view,'Allen brain atlas ontology' ) ## brain regions names

        self.right_panels.addWidget(self.trees_view)

        self.setRightPanelStrechFactor()
        self.frame_splitter.setStretchFactor(0, 8)
        self.frame_splitter.setStretchFactor(1, 1)
        self.frame_splitter.setCollapsible(0, False)
        self.frame_splitter.setCollapsible(1, True)

        # add frame_splitter to main layout
        self.MainLayout.addWidget(self.frame_splitter)

    def SynAtlas(self):
        print('syn in action')

    def setRightPanelStrechFactor(self):
        self.right_panels.setStretchFactor(0, 4)
        self.right_panels.setStretchFactor(1, 5)
        self.right_panels.setStretchFactor(2,0)
        self.right_panels.setCollapsible(0, True)
        self.right_panels.setCollapsible(1, True)
        self.right_panels.setCollapsible(2, True)
        self.right_panels.setSizes([360, 360, 200])

    def add_menubar(self):
        """
        set up manu bar for the main window
        """
        self.mbar = pg.QtGui.QMenuBar()
        self.MainLayout.setMenuBar(self.mbar)
        self.OptionMenu = self.mbar.addMenu('&Tools')
        self.OptionAction1 = pg.QtGui.QAction("&Pickle current images")
        self.OptionAction1.triggered.connect(self.pickeImage_clicked)
        self.OptionMenu.addAction(self.OptionAction1)
        self.OptionAction2 = pg.QtGui.QAction("&Batch downloading")
        self.OptionAction2.triggered.connect(self.batchDownloading_clicked)
        self.OptionMenu.addAction(self.OptionAction2)

        self.HelpMenu = self.mbar.addMenu('&Help')
        self.HelpAction = pg.QtGui.QAction("&LICENSE")
        self.HelpAction.setStatusTip('BSD-3')
        self.HelpAction.triggered.connect(self.License_clicked)
        self.HelpMenu.addAction(self.HelpAction)

    def imageHoverEvent(self, event):
        """Show the position, pixel, and value under the mouse cursor.
        """
        if event.isExit():
            self.Expression_view.view.setTitle("")
            self.ISH_view.view.setTitle("")
            return

        # data = self.ISH_view.image
        if type(self.Expression_view.image)!=type(None):
            if self.ISH_view.view.sceneBoundingRect().contains(event.pos().x(), event.pos().y()):
                ppos = self.ISH_view.view.mapToView(event.pos())
                i, j = int(ppos.x()), int(ppos.y())
                if i >=0 and j >=0 and i < self.Expression_view.image.shape[1] and j <self.Expression_view.image.shape[0]:
                    fishVal = self.ISH_view.image[j,i, :]
                    self.ISH_view.view.setTitle(f"pos:({i},{j}), rgb: {fishVal}")
                    if len(self.Expression_view.image.shape)==3:
                        expressionVal = self.Expression_view.image[j, i, :]
                        pixelIntensity = rgb2intensity(np.reshape(expressionVal,(1,1,3)), jetcm)[0][0]
                        rgb_txt = f'RGB: {expressionVal}.'

                    else:
                        rgb_txt = ''
                        pixelIntensity = self.Expression_view.image[j, i]

                    if self.atlas_ImageID!=None:
                        expTitle = f"pos:({i},{j}),"+rgb_txt +f" Expresssion intensity:{pixelIntensity:.2f}. AtlasImage ID: {self.atlas_ImageID}"
                    else:
                        expTitle =  f"pos:({i},{j}),"+rgb_txt +f" Expresssion intensity:{pixelIntensity:.2f}. AtlasImage ID: {self.atlas_ImageID}"

                    self.Expression_view.view.setTitle(expTitle)

                    self.ISH_view_vLine.setPos(ppos.x())
                    self.ISH_view_hLine.setPos(ppos.y())
                    self.ISH_view_vLine.setZValue(1000)
                    self.ISH_view_hLine.setZValue(1000)

                    self.Expression_view_vLine.setPos(ppos.x())
                    self.Expression_view_hLine.setPos(ppos.y())
                    self.Expression_view_vLine.setZValue(1000)
                    self.Expression_view_hLine.setZValue(1000)

                else:
                    self.Expression_view.view.setTitle("")
                    self.ISH_view.view.setTitle("")

    def License_clicked(self):
        with open("Data/LICENSE.txt") as f:
            BSD_3 = f.readlines()
        BSD_3b = ""
        for l in BSD_3:
            BSD_3b = BSD_3b+l
        showdialog("This program is under BSD-3 license.\nCopyright (c) 2020-2021, ZeitgeberH@github. All rights reserved.", BSD_3b,False)

    def batchDownloading_clicked(self):
        WindowTemplate, TemplateBaseClass = pg.Qt.loadUiType('Data/uis/batchDownloadOptionsDialog.ui')
        self.bd_dialog  = TemplateBaseClass()
        self.bd_dialog_form = WindowTemplate()
        self.bd_dialog_form.setupUi(self.bd_dialog)
        self.bd_dialog_form.fileButton.clicked.connect(self.batchGenesFile)
        self.bd_dialog_form.buttonBox.clicked.connect(self.buttonBoxResponse)
        self.bd_dialog.exec_()

    def getBatchDownloadInfo(self,showDialog=False):
        fullResolution =self.bd_dialog.children()[2].children()[0].children()[1].isChecked()
        Downsampled =self.bd_dialog.children()[2].children()[0].children()[2].isChecked()
        downSampledFactor = self.bd_dialog.children()[2].children()[0].children()[0].children()[4].value()
        QFactor = self.bd_dialog.children()[2].children()[0].children()[0].children()[5].value()

        for j in range(3):
            if self.bd_dialog.children()[2].children()[1].children()[j].isChecked():
                planeOption = j
                break
        userDefinedGenesFiles = self.bd_dialog_form.filePath.text()
        with open(userDefinedGenesFiles,'r') as file:
            geneNames = file.read().splitlines()
        geneNames = [g for g in geneNames if len(g)>0]
        self.BatchDownloading_FileInfo(geneNames, planeOption, fullResolution, downSampledFactor, QFactor,showDialog)

    def buttonBoxResponse(self,evt):
        '''
        if 'OK' button pressed, parse values from the UI

        Parameters
        ----------
        evt : TYPE: buttonBox responses
            DESCRIPTION.

        Returns
        -------
        None.

        '''
        if evt.text() == 'OK':
            if self.bd_dialog_form.filePath.text()!='':
                self.getBatchDownloadInfo(True)
            else:
                showdialog("Press 'Locate file' button to specify the file with names of genes you would like to know!")

    def BatchDownloading_FileInfo(self, geneNames, planeOption, fullResolution, downSampledFactor, QFactor, showDialog=False):
        p = ['all','coronal','sagittal']
        if fullResolution:
            ds = f'\nDownloading {p[planeOption]} planes section images for following genes at full resolution:'
        else:
            ds = f'\nDownloading {p[planeOption]} section images at downsampling factor:{downSampledFactor}, quality:{QFactor} for following genes at reduced of {planeOption}: '
        print(ds)
        print(f'{geneNames}')

        geneDict = {}
        inforStr = 'Gene    Experiments\n'
        for g in geneNames:
            geneExpID_all, geneExpData_all = self.geapi.get_gene_experiments2(g)
            geneExpID = []
            geneExpData = []
            if geneExpID_all!=None:
                if planeOption!=0:
                    for gid, ged in zip(geneExpID_all,geneExpData_all):
                        if ged['plane_of_section_id']== planeOption: ## filter specified image planes!
                            geneExpID.append(gid)
                            geneExpData.append(ged)
                else:
                    geneExpID = geneExpID_all
                    geneExpData = geneExpData_all

                inforStr = inforStr+ f'{g:8}  {len(geneExpID)}\n'
                geneDict[g]={'geneExpID': geneExpID, 'geneExpData':geneExpData}

        returnValue = showdialog(ds, details= inforStr)

        if returnValue == QtWidgets.QMessageBox.Ok:
            self.BatchDownload(geneDict, fullResolution, downSampledFactor, QFactor)

    def BatchDownload(self,geneDict, imageQuality, downSampleFactor, qualityFactor):
        Expression=False
        with pg.ProgressDialog("Downloading FISH images for gene..." , maximum=len(geneDict),\
                               busyCursor=True, nested=True) as dlg_gene:
            for gene in geneDict:
                expIDs = geneDict[gene]['geneExpID']
                with pg.ProgressDialog("Experiments..." , maximum=len(expIDs),\
                                       busyCursor=True, nested=True) as dlg_exp:
                    for expID in expIDs:
                        cache = check_gene_cached_images(self.geapi.fish_images_cache, gene, expID,Expression, imageQuality, downSampleFactor, qualityFactor)
                        #### check cache and download if not availabe
                        if not cache:
                            imageIDs, sectionNumbers, expMetaData = self.geapi.get_gene_experiments_imageList(expID) ## query imageIDs
                            with pg.ProgressDialog("... ISH & Expression" , maximum=2,\
                                                   busyCursor=True, nested=True) as dlg_eh:

                                for Expression in [True, False]: ## expression images
                                    with pg.ProgressDialog("...... Images..." , maximum=len(imageIDs),\
                                                           busyCursor=True, nested=True) as dlg1:
                                        for j, sn in zip(imageIDs, sectionNumbers):
                                            self.geapi.download_gene_section_imageData(gene, expID, j, sn, Expression, imageQuality, downSampleFactor,qualityFactor)
                                            dlg1 += 1
                                            if dlg1.wasCanceled():
                                                print("Canceled stage")
                                                break
                                    dlg_eh += 1
                                    if dlg_eh.wasCanceled():
                                        print("Canceled images stage")
                                        break
                        dlg_exp += 1
                        if dlg_exp.wasCanceled():
                            print("Canceled experiment stage")
                            break

                dlg_gene += 1
                if dlg_gene.wasCanceled():
                    print("Canceled gene stage")
                    break

    def batchGenesFile(self):
        fileName = getfiles()
        self.bd_dialog_form.filePath.setText(fileName)
        # self.getBatchDownloadInfo(showDialog=True)

    def add_toolbar(self):
        self.toolbar = QtWidgets.QToolBar("Main toolbar")
        self.addToolBar(2, self.toolbar) # https://doc.qt.io/qt-5/qt.html#ToolBarArea-enum
        self.toolbar.addSeparator()
        self.toolbar.addSeparator()
        # self.tooglesettingsAction = pg.QtGui.QAction(pg.QtGui.QIcon("Data/icons/lsm.png"), "Open LSM images")
        # self.tooglesettingsAction.triggered.connect(self.openLSM_clicked)
        # self.toolbar.addAction(self.tooglesettingsAction)

    def sycImageToRegion(self,expID, regionName):
        ##TODO
        try:
            imgData_meta = self.geapi.getExpTargetXY(expID, regionName)
            print('current region:  '+ regionName)

            imgPath = os.path.join(self.currentImageCacheDir, str(imgData_meta['section_number']) + '_'+str(imgData_meta['section_image_id']))
            imgPath_E = imgPath+'_E.jpg'
            imgPath_H = imgPath+'_H.jpg'
            self.updateFishImage(self.ISH_view, imgPath_H)
            self.updateFishImage(self.Expression_view, imgPath_E, True)
            row = self.currentSearchISH_imageLst.index(imgPath_H)
            self.experimentImage_Tab.setCurrentCell(row, 0)
            self.experimentImage_Tab.setFocus()
            return True

        except:
            showdialog('Syncronization failed!')
            return False

    def event_parameters_stateChange(self, params, changes):
        '''
        Responding to user's input in the parameter tree

        Parameters
        ----------
        params : TYPE
            DESCRIPTION.
        changes : TYPE
            DESCRIPTION.

        Returns
        -------
        None.

        '''
        pv = self.globalPars.getValues()
        for param, change, data in changes:
            childName = param.name()
            if childName == 'Gene':
                self.experiments_view.setCurrentIndex(0)

            if childName=='region':
                if self.currentSearchDict!=None:
                    self.currentSearchDict['region'] = data
            if (childName == 'Sync') and (self.currentExpID!=None):
                regionName = pv['Gene & Brain region & species'][1]['Region'][0]
                if regionName!='':
                    print(f'Sync to region {regionName}')
                else:
                    showdialog("Specifiy a rgion to syn!")
            if (childName == 'Syn to atlas after double clicking'):
                if not data:
                    self.clearSVGItems() ## turn off atlas overlay
            if childName=='Expression':
                if data=='color mask':
                    self.updateFishImage(self.Expression_view,self.Expression_mask_file, True, False)
                else:
                    self.updateExpressionIntensity()
            if childName == 'atlas section offset':
                # pdb.set_trace()
                if self.atlasID!=None:
                    idx = self.svgItem_atlasE_index
                    atlas_ofs_idx = pv['image-to-Atlas'][1]['atlas section offset'][0] + idx
                    if atlas_ofs_idx>-1 and atlas_ofs_idx<len(self.atlas_sectionID_list):
                        atlas_imageID = self.atlas_sectionID_list[atlas_ofs_idx]
                        atlas_fn  = self.geapi.getAtlasBoundaryImage(atlas_imageID, self.atlasID, GraphicGroupLabel_id=28)
                        xP, yP = self.svgItem_atlasE.pos().x(),self.svgItem_atlasE.pos().y()
                        self.overlayAtlasSVG(xP+self.atlas_x,yP+self.atlas_y,self.atlas_x,self.atlas_y, atlas_fn)
                        self.atlas_ImageID = atlas_imageID
                else:
                    print('Double clicking to syn an atlas image first!')

            if childName == 'x offset':
                # pdb.set_trace()
                x, y = self.svgItem_atlasE_pos
                y_offset = pv['image-to-Atlas'][1]['y offset'][0]
                self.svgItem_atlasE.setPos(x+data, y+y_offset)
                # self.svgItem_atlasH.setPos(x+data, y+y_offset)

            if childName == 'y offset':
                # pdb.set_trace()
                x_offset = pv['image-to-Atlas'][1]['x offset'][0]
                x, y = self.svgItem_atlasE_pos
                self.svgItem_atlasE.setPos(x+x_offset, y+data)
                # self.svgItem_atlasH.setPos(x+x_offset, y+data)

            if childName == 'Expression mask opacity':
                opc = pv['image-to-Atlas'][1][childName][0]
                self.Expression_view.view.setOpacity(opc)

        if (len(changes) == 1) and (childName == 'Sync') and (self.currentExpID!=None):
            if self.geapi!=[]:
                regionName = pv['Gene & Brain region & species'][1]['Region'][0]
                self.sycImageToRegion(self.currentExpID, regionName.upper())
            else:
                showdialog("Load an experiment first!")
            return

        #
        search_pars = {}
        search_pars['gene'] = pv['Gene & Brain region & species'][1]['Gene'][0]
        if search_pars['gene']=='':
            return
        search_pars['region'] = pv['Gene & Brain region & species'][1]['Region'][0]
        search_pars['species'] = pv['Gene & Brain region & species'][1]['Species'][0]
        search_pars['plane'] = pv['View options'][1]['Plane'][0]
        search_pars['Syn to atlas after double clicking'] = pv['image-to-Atlas'][1]['Syn to atlas after double clicking'][0]
        search_pars['maximal resolution'] = pv['View options'][1]['Maximal resolution'][0]
        search_pars['downsample'] = pv['View options'][1]['Image quality at lower resolution'][1]['Downsample'][0]
        search_pars['quality'] = pv['View options'][1]['Image quality at lower resolution'][1]['Downsample'][0]
        self.updateGeneOfInterest(search_pars)

    def resetMergeImage(self):
        for item in self.lsm_bottomRight_view.view.items:
            if hasattr(item, 'Name'):
                if item.Name[:3]!='br_':
                    self.lsm_bottomRight_view.view.removeItem(item)
        self.lsm_bottomRight_view.setImage(self.mergeImage, axes={'x':0,'y':1,'c':2,'t':None},autoRange=False,autoLevels=False)

    def overlay_segmentation(self, bds, imageViewHandle, noMask=False):
        '''
        Overlay segmented ROIs to given ImageView instance

        Parameters
        ----------
        bds : TYPE: ndarray
            contour boundaries of segemented results
        imageViewHandle : pyqtgraph.imageView instance
            imageview window to draw

        Returns
        -------
        None.
        '''
        for item in imageViewHandle.view.items:
            if hasattr(item, 'Name'):
                if item.Name[:3]=='sgb':
                    imageViewHandle.view.removeItem(item)
        if noMask: ## no mask
            return
        if len(bds.shape)==2:
            imItem = pg.ImageItem(bds)
            imItem.Name = 'sgb'
            imageViewHandle.view.addItem(imItem)
            imItem.setZValue(100)
        else:
            allChans = set([0,1,2])
            for nMask in range(bds.shape[2]):
                mask_ =  np.empty_like(bds)
                mask_[:,:,nMask] = bds[:,:,nMask]
                otherChans = list(allChans.difference([nMask]))
                for j in otherChans:
                    mask0 =  np.empty_like(bds[:,:,0])
                    mask0[bds[:,:,nMask]>0] = 80 ## set other channel's pixels intensity. so we can get a color image
                    mask_[:,:,j] = mask0
                imItem = pg.ImageItem(mask_)
                imageViewHandle.view.addItem(imItem)
                imItem.Name = 'sgb'+str(nMask)
                imItem.setZValue(100)

    def initGeneAPI(self):
        try:
            self.geapi = GeneRenderAPI()
            self.geapi.DOWNSAMPLE_FACTOR = 4 ## downsampling factor if not Full image, range(0 to 10), higher the lower quality
            self.geapi.QUALITY_FACTOR = 50 ## quality factor if not Full image, range(0 to 100), higher the better quality
            return True
        except:
            self.geapi = []
            print('Error to initiate gene query API. Try again')
            return False

    def loadSectionImageInfo(self,row, col):
        '''
        load section image metainformation

        Parameters
        ----------
        row : TYPE
            DESCRIPTION.
        col : TYPE
            DESCRIPTION.

        Returns
        -------

        '''
        expID = self.experimentsInfo_Tab.item(row, 0).value
        expPlane = self.experimentsInfo_Tab.item(row, 1).value
        expGene = self.experimentsInfo_Tab.item(row, 2).value

        if self.geapi==[]:
            r=self.initGeneAPI()
            if not r:
                return
        ## query and return section image list from this experiment
        ##TODO
        ## CHECK cached image file first.If no, then request from ABA
        imageIDs, secNums, imageMeta = self.geapi.get_gene_experiments_imageList(expID)
        expDF = pd.DataFrame()
        expDF.loc[:,'Section'] = secNums
        expDF.loc[:,'imageID'] = imageIDs
        expDF['Experiment'] = expID
        expDF['Gene'] = expGene
        expDF['Plane'] = expPlane
        expDF['Height'] = [e['image_height'] for e in imageMeta]
        expDF['Width'] = [e['image_width'] for e in imageMeta]
        expDF['Resolution'] = [e['resolution'] for e in imageMeta]
        expDF['Failed'] = [e['failed'] for e in imageMeta]
        expDF = expDF.set_index('Section')
        if expPlane == 'coronal':
            expDF.sort_index(ascending=False, inplace=True)
        data = np.array(expDF.to_records(),
                              dtype=[('Section',object), ('imageID', object), ('Experiment', object), ('Gene', object),\
                                     ('Plane', object),('Height', object),('Width', object),('Resolution', object),\
                                        ('Failed', object) ])
        self.experimentImage_Tab.setData(data)
        self.experiments_view.setCurrentIndex(1)
        self.currentExpID = expID
        # self.currentSearchISH_imageLst = imageFilesList_H
        # self.currentSearchExpression_imageLst = imageFilesList_E
        self.loadExp_image(expGene, expID, expPlane,0)

    def render_sectionImage(self,row, col, previousRow, previousColumn):
        '''
        render section images in imageItems panels

        Parameters
        ----------
        row : TYPE
            DESCRIPTION.
        col : TYPE
            DESCRIPTION.

        Returns
        -------
        None.

        '''
        if self.currentSearchISH_imageLst !=None:
            if self.currentISH_image != self.currentSearchISH_imageLst[row]:
                self.clearSVGItems()
                self.updateFishImage(self.ISH_view, self.currentSearchISH_imageLst[row])
                self.Expression_mask_file = self.currentSearchExpression_imageLst[row]
                self.Expression_intensityVal = None
                self.updateFishImage(self.Expression_view, self.currentSearchExpression_imageLst[row], True)
                self.Expression_view.view.setXLink(self.ISH_view.view)
                self.Expression_view.view.setYLink(self.ISH_view.view)
                self.currentISH_image = self.currentSearchISH_imageLst[row]
                self.statusBar.showMessage(self.currentISH_image , 5000)

    def loadExp_image(self, gene, expID,plane, idx=0):
        '''
        Load FISH image list for current gene and experiment ID with image index at index into table
        If not in cache, download it from ABI.

        Parameters
        ----------
        gene : TYPE, str
            Gene symbol.
        expID : TYPE, Int
            Experiment ID
        idx : TYPE, optional
            Section index to load for current experiment. The default is 0.
        fullImage : TYPE, bool
            DESCRIPTION. whether to download the full resolution images. Would take times if set true (and possible connection failure)

        Returns
        -------
        None.

        '''
        imageQuality = self.currentSearchDict['maximal resolution']
        downSampleFactor =self.currentSearchDict['downsample']
        qualityFactor = self.currentSearchDict['quality']
        Expression=False
        ##TODO
        ## check by image IDs. (for now, it just check if there is ANY images)
        cache = check_gene_cached_images(self.geapi.fish_images_cache, gene, expID,Expression, imageQuality, downSampleFactor, qualityFactor)
        #### check cache and download if not availabe
        if not cache:
            self.statusBar.showMessage(f'Downloading data for gene {gene}. Please wait...', 3000)
            imageIDs, sectionNumbers, expMetaData = self.geapi.get_gene_experiments_imageList(expID) ## query imageIDs
            with pg.ProgressDialog("Downloading FISH images..." , maximum=2,\
                                   busyCursor=True, nested=True) as dlg0:

                for Expression in [True, False]: ## expression images
                    with pg.ProgressDialog("Downloading ..." , maximum=len(imageIDs),\
                                           busyCursor=True, nested=True) as dlg1:
                        for j, sn in tqdm.tqdm(zip(imageIDs, sectionNumbers), total=len(imageIDs)):
                            self.geapi.download_gene_section_imageData(gene, expID, j, sn, Expression, imageQuality, downSampleFactor,qualityFactor)
                            dlg1 += 1
                            if dlg1.wasCanceled():
                                print("Canceled stage")
                                break
                        dlg0 += 1
                        if dlg0.wasCanceled():
                            print("Canceled stage")
                            break
                else:
                    cache=True

        if cache: ## load cached data
            imageFilesList_H = self.geapi.get_gene_experiments_imageData(gene, expID,True, False, imageQuality,downSampleFactor, qualityFactor)## ISH IMAGES
            # self.updateFishImage(self.ISH_view, imageFilesList_H[idx])
            imageFilesList_E = self.geapi.get_gene_experiments_imageData(gene, expID,True, True,imageQuality,downSampleFactor, qualityFactor)## EXPRESSION MASKS
            # self.updateFishImage(self.Expression_view, imageFilesList_E[idx])
            if plane == 'coronal':
                imageFilesList_H  = imageFilesList_H [::-1]
                imageFilesList_E = imageFilesList_E[::-1]
            self.currentSearchISH_imageLst = imageFilesList_H
            self.currentSearchExpression_imageLst = imageFilesList_E
            self.currentImageCacheDir = str(Path(imageFilesList_E[0]).parent)
            region = self.globalPars.getValues()['Gene & Brain region & species'][1]['Region'][0]
            self.sycImageToRegion(self.currentExpID, region.upper())
        else:
            print('Glitches in cached FISH images...')

    def clearSVGItems(self):
        self.atlas_ImageID = None
        for item in self.Expression_view.view.items:
            if type(item) == QtSvg.QGraphicsSvgItem:
                self.Expression_view.view.removeItem(item)
                item.setParent(None)
                del item
        for item in self.ISH_view.view.items:
            if type(item) == QtSvg.QGraphicsSvgItem:
                self.ISH_view.view.removeItem(item)
                item.setParent(None)
                del item

    def otherInitStuff(self):
        self.initGeneAPI()
        self.initABAtlasTable()
        self.currentSearchDict = None
        self.currentSearchISH_imageLst = None
        self.currentSearchExpression_imageLst = None
        self.currentExpID = None
        self.currentISH_image = None
        self.defaultBrowser = None
        self.atlas_sectionID_list = None
        self.atlas_x = None
        self.atlas_y = None
        self.atlasID = None
        self.atlas_ImageID = None
        self.Expression_mask = None
        self.Expression_intensity = None
        ## LSM file related
        # self.currentLSMfile = None
        # self.currentLSMimage = None
        # self.segs_merge = {'R':None,'G':None,'B':None} ## store segmented value right before postprocessing steps
        # self.merge_r = True
        # self.merge_g = True
        # self.merge_b = True
        # self.DNNModel = None ## deal with tf.function retracing?

    def initABAtlasTable(self):
        ## load Allen Brain isntitute atlas acronym
        df = self.geapi.atlas_df[['id','acronym','name']]
        self.abAtlas_view.setData(df.to_records())

    def getSynedAtlas(self,atlasID, imageID, x,y, groupID=28):
        imgSynDict = self.geapi.synImageToAtlas(atlasID, imageID, (x, y)) ## get  dictionary of syned atlas Image information
        #imgSynDict['section_image_id'],imgSynDict['x'],imgSynDict['y']
        if len(imgSynDict)>0:
            atlas_imageID, xa,ya = imgSynDict['section_image_id'], imgSynDict['x'],imgSynDict['y']
            print(f'Syn info: section imgeID {imageID}, pos:({x},{y}). Atlas imageID {atlas_imageID}, pos: ({xa},{ya})')
            fileName  = self.geapi.getAtlasBoundaryImage(atlas_imageID, atlasID, GraphicGroupLabel_id=groupID)
            return fileName, atlas_imageID,xa, ya
        return None

    def overlayAtlasSVG(self, x, y, xa,ya, atlas_fileName):
        '''
        Parameters
        ----------
        x : TYPE, float
            current section image point's x
        y : TYPE, float
            current section image point's x
        xa : TYPE, float
            atlas image synced X
        ya atlas_fileName : TYPE
            atlas image synced Y
        Returns
        -------
        None.

        '''

        self.clearSVGItems()
        self.svgItem_atlasE = QtSvg.QGraphicsSvgItem(atlas_fileName)

        self.svgItem_atlasE.setOpacity(1);#0.08)
        self.svgItem_atlasE.scale(1, 1)
        self.svgItem_atlasE.setPos(x-xa, y-ya)
        self.svgItem_atlasE_pos = (x-xa, y-ya) ## for reference
        self.Expression_view.view.addItem(self.svgItem_atlasE)
        self.svgItem_atlasE.setZValue(-100)

    def ISH_view_mouseClicked(self, event):
        if type(self.ISH_view.image)!=type(None):
            if self.currentSearchDict['Syn to atlas after double clicking']: ## only syn when True
                if event.button()==1 and event.double(): ## if left double clicking~
                    if self.ISH_view.view.sceneBoundingRect().contains(event.pos().x(), event.pos().y()):
                            if self.currentSearchDict['species']=='Mouse':
                                row = self.experimentImage_Tab.currentRow()
                                imageID = self.experimentImage_Tab.item(row, 1).value ## get current section image ID
                                planeName = self.experimentImage_Tab.item(row, 4).value
                                if planeName == 'coronal':
                                    atlasID = 1
                                else:
                                    atlasID = 2
                                x = self.ISH_view_vLine.pos().x() ## position value from event is not consistent with hover event!!!
                                y = self.ISH_view_hLine.pos().y()

                                atlas_fn,atlas_imageID, xa,ya = self.getSynedAtlas(atlasID, imageID, x, y)
                                self.overlayAtlasSVG(x,y,xa,ya, atlas_fn)
                                ##TODO
                                if self.atlas_sectionID_list==None:
                                    self.atlas_sectionID_list, _,_ = self.geapi.getAtlasImagesAnnoationList(atlasID)
                                self.svgItem_atlasE_index  = self.atlas_sectionID_list.index(atlas_imageID) ## store it for manul correction
                                self.atlas_x = xa ## same above
                                self.atlas_y = ya
                                self.atlasID = atlasID
                                self.atlas_ImageID = atlas_imageID
                            else:
                                print('Not supported yet for human data')

    def updateExperimentsTable(self, expDF):
        data = np.array(expDF.to_records(),
                              dtype=[('Experiment', object),('Plane', object), ('Gene symbol', object),\
                                     ])
        self.experimentsInfo_Tab.setData(data)

    def updateGeneOfInterest(self, parDict):
        if self.geapi==[]:
            r=self.initGeneAPI()
            if not r:
                return
        if parDict['species']=='Mouse':
            parDict['gene'] = parDict['gene'].capitalize()
        else:
            parDict['gene'] = parDict['gene'].upper()
        self.updateGenessTable(parDict['gene'],  parDict['species'])
        try:
            geneExpID, geneExpData = self.geapi.get_gene_experiments2(parDict['gene'])
        except:
            print('No respondes from Server!')
            return
        if geneExpID != None:
            expDF = pd.DataFrame()
            expDF.loc[:,'Experiment'] = geneExpID
            expDF.loc[:,'Plane'] = [self.geapi.sectionPlaneDict[d['plane_of_section_id']] for d in geneExpData]
            expDF['Gene Symbol'] = parDict['gene']
            # expDF.loc[:,'Expression Summary'] = [] ## TODO
            expDF = expDF.set_index('Experiment')
            if parDict['plane']!='All':
                expDF=expDF[expDF['Plane']==parDict['plane'].lower()]
            self.updateExperimentsTable(expDF)
            self.statusBar.showMessage(f'Data for gene {parDict["gene"]} is listed!', 10000)
            self.currentSearchDict = parDict ## updating state variables!
            self.geapi.DOWNSAMPLE_FACTOR = parDict['downsample']
            self.geapi.QUALITY_FACTOR = parDict['quality']
        else:
            self.statusBar.showMessage(f'Data for gene {parDict["gene"]} is not found', 10000)

    def updateExpressionIntensity(self):
        # self.Expression_mask
        with pg.BusyCursor():
            if self.Expression_intensityVal is None:
                if self.Expression_view.image is None:
                    img = np.array(Image.open(self.Expression_mask_file))
                else:
                    img = self.Expression_view.image
                self.Expression_intensityVal = rgb2intensity(img)
                self.Expression_view.setImage(self.Expression_intensityVal, pos = [0,0], axes={'x':1,'y':0,'c':None,'t':None},autoRange=False)
            else:
               self.Expression_view.setImage(self.Expression_intensityVal, pos = [0,0], axes={'x':1,'y':0,'c':None,'t':None},autoRange=False)

    def updateFishImage(self, imageItemHandle, imgFileName, Expression=False, autoRange=True):
        '''set image for imageView
        '''
        # imgFileName = base_dir / 'FishImagesCache' / 'Stac-69887313' / 'DQ' /'65_69866246_H.jpg'
        # img = np.transpose(Image.open(imgFileName),axes=(1,0,2))
        img = np.array(Image.open(imgFileName))
        imageItemHandle.clear()
        # w,h,_ = img.shape
        # imageItemHandle.setImage(img, pos = [h/2,w/2])
        imageItemHandle.setImage(img, pos = [0,0], axes={'x':1,'y':0,'c':2,'t':None},autoRange=autoRange)
        if Expression:
            imageItemHandle.view.setOpacity(0.85)
        # imageItemHandle.setOpts(update=True, opacity=0.2)

    def findGenesInTable(self,geneSymbol):
        nrows = self.geneInfo_Tab.rowCount()
        for r in range(nrows):
            if geneSymbol==self.geneInfo_Tab.item(r,0).value:
                return True
        else:
            return False

    def updateGenessTable(self, geneSymbol, species):
        ##TODO
        #### CHECK if gene is in cache
        if not self.findGenesInTable(geneSymbol):
            geneInfo = queryGene_ncbi(geneSymbol, species)
            if len(geneInfo) > 0:
                expDF = pd.DataFrame([geneInfo])
                expDF['symbol']= geneSymbol
                expDF['genecard url'] = 'genecards.org/cgi-bin/carddisp.pl?gene='+geneSymbol
                expDF = expDF[['symbol','description','type',\
                       'entrez_id','chromosomes', 'number of transcripts','genecard url']]
                data = np.array(expDF.to_records(index=False),
                                      dtype=[('symbol', object),('description', object), ('type', object),\
                                             ('entrez_id', object),('chromosomes', object),('number of transcripts', object),('genecard url', object)
                                             ])
                self.geneInfo_Tab.appendData(data)

    def querySelectedGeneInfo(self, row, col):
        if col== self.geneInfo_Tab.columnCount()-1:  ## last col
            url =  self.geneInfo_Tab.item(row, col).value
            self.openGeneCardPage(url)

    def openGeneCardPage(self, url):
        if self.defaultBrowser==None: ##  if default browser not set, try get default
            try:
                from gr_allenFISH_utils import getDefaultBrowerPath
                browserPath = getDefaultBrowerPath().replace('\\','/') +' %s'
                self.defaultBrowser = browserPath
                webbrowser.get(using = self.defaultBrowser).open_new(url)
            except:
                webbrowser.open(url,new=2) ## fallback to IE if failed
        elif self.defaultBrowser == -1: ## failed already
            webbrowser.open(url,new=2) ## fallback to IE
        else: ## properly set
            webbrowser.get(using = self.defaultBrowser).open_new(url)

    def pickeImage_clicked(self):
        if self.currentExpID is not None:
            import pickle
            with pg.BusyCursor():
                fname = self.currentSearchDict['gene']+'_'+self.currentISH_image+'_ExpID-'+str(self.currentExpID)
                fishVal = self.ISH_view.image
                expressionVal = self.Expression_view.image
                with open(fname+'_expression.pickle','wb') as handle:
                    pickle.dump(expressionVal, handle, protocol=pickle.HIGHEST_PROTOCOL)

                with open(fname+'_fish.pickle','wb') as handle:
                    pickle.dump(fishVal,handle, protocol=pickle.HIGHEST_PROTOCOL)
        else:
            print('No image loaded yet!')

    def closeFISH_ViewItems(self):
        for view in [self.ISH_view,self.Expression_view]:
            view.clear()
            view.close()

    def closeEvent(self, event):
        reply = QMessageBox.question(self, 'Closing Window', 'Are you sure you want to close the window?',\
                                     QMessageBox.Yes | QMessageBox.No, QMessageBox.No)
        if reply == QMessageBox.Yes:
            try: ## release resources
                self.closeFISH_ViewItems()
            except:
                print('clean up failed!')
            event.accept()
        else:
            event.ignore()
def main():
    app = QtWidgets.QApplication(sys.argv)
    main = MainWindow(app)
    main.show()
    sys.exit(app.exec_())

if __name__ == '__main__':
    main()

max_stars_repo_name: ZeitgeberH/FISH-VIEWER

max_issues_repo_name: ZeitgeberH/FISH-VIEWER

max_forks_repo_name: ZeitgeberH/FISH-VIEWER


text:
import math
import sys
import numpy as np

class TransferProbability:
    def __init__(self, network):
        self.nodes = network.nodes
        self.edges = network.edges
        self.trajectories = network.trajectories
        self.derive()

    def derive(self):
        """
        Derive matrix P ,matrix Q ,matrix S and column vector V of each node.

        """
        for node in self.nodes:
            p = self.create_transition_matrix(node)
            q, s = self.reorganize(p, node)
            node.vector = self.cal_vector(node, p, q)

    def create_transition_matrix(self, d):
        """
        Construct the transition matrix P by function transition_probability.

        :param Point d: the destination node
        :return: matrix P

        """
        nodes_len = len(self.nodes)
        p = np.zeros((nodes_len, nodes_len), dtype=float)
        p_row, p_col = p.shape[0], p.shape[1]
        for row in range(p_row):
            for col in range(p_col):
                p[row][col] = self.transition_probability(d, self.nodes[row],
                                                          self.nodes[col])

        return p

    def transition_probability(self, d, nodei, nodej):
        """
        Get the transition probability of moving from nodei to nodej
        through the state of nodei and the subscripts of both nodes.

        :param Point d: the destination node
        :param Point nodei: the starting node of transition
        :param Point nodej: the ending node of transition
        :return: the transition probability

        """
        absorbing_state = [-1 for j in range(len(self.edges))]
        nodei_index = self.nodes.index(nodei)
        nodej_index = self.nodes.index(nodej)

        if (nodei == d or self.edges[nodei_index] == absorbing_state) \
                and nodei_index == nodej_index:
            return 1
        elif (not (nodei == d or self.edges[nodei_index] == absorbing_state)) \
                and nodei_index != nodej_index:
            return self.prd(d, nodei, nodej)
        else:
            return 0

    def prd(self, d, nodei, nodej):
        """
        Get the turning probability of moving from nodei to nodej
        through the ratio of adding func values of all the trajectories
        on (nodei,nodej) and all the trajectories starting from nodei.

        :param Point d: the destination node
        :param Point nodei: the starting node of transition
        :param Point nodej: the ending node of transition
        :return: the turning probability

        """
        sum_ij, sum_i = 0, 0
        nodei_index = self.nodes.index(nodei)
        nodej_index = self.nodes.index(nodej)

        # add func values of all the trajectories on (nodei,nodej)
        if self.edges[nodei_index][nodej_index] != -1:
            for t in self.edges[nodei_index][nodej_index]:
                sum_ij += self.func(t, d, nodei)

        # add func values of all the trajectories starting from nodei
        for col in range(len(self.edges)):
            if self.edges[nodei_index][col] != -1:
                for t in self.edges[nodei_index][col]:
                    sum_i += self.func(t, d, nodei)

        if sum_i == 0:
            return 0
        return sum_ij / sum_i

    def func(self, traj, d, nodei):
        """
        Estimate the likelihood that a trajectory traj might
        suggest a correct route to d.

        :param Point.trajectory_id traj: the trajectory
        :param Point d: the destination node
        :param Point nodei: the starting node
        :return: the likelihood

        """
        dists = sys.maxsize
        flag = 0
        traj_value = self.trajectories[traj]
        nodei_index_in_traj = \
            self.trajectories[traj].index(self.nodes.index(nodei))

        # trajectory traj passes node d ,dists = 0
        if self.nodes.index(d) in traj_value[nodei_index_in_traj::]:
            dists = 0
            flag = 1
        # trajectory traj only has one node rather than edge
        elif len(traj_value[nodei_index_in_traj::]) == 1:
            dists = math.pow(
                ((d.latitude - self.nodes[traj_value[0]].latitude) ** 2
                 + (d.longitude - self.nodes[traj_value[0]].longitude) ** 2),
                0.5)
            flag = 1
        # trajectory traj has one edge at least
        elif len(traj_value[nodei_index_in_traj::]) >= 2:
            for index in range(nodei_index_in_traj, len(traj_value) - 1):
                new_dist = self.get_dist(d, self.nodes[traj_value[index]],
                                         self.nodes[traj_value[index + 1]])
                if new_dist < dists:
                    dists = new_dist
                    flag = 1

        if flag == 0:
            return 0
        return math.exp(-dists)

    def get_dist(self, d, point1, point2):
        """
        Get the shortest Euclidean/network distance between d and
        the segment from point1 to point2.

        :param Point d: the point outside the segment
        :param Point point1: the endpoint of the segment
        :param Point point2: another endpoint of the segment
        :return: the distance

        """
        d_x = d.latitude
        d_y = d.longitude
        point1_x = point1.latitude
        point1_y = point1.longitude
        point2_x = point2.latitude
        point2_y = point2.longitude
        cross = (point2_x - point1_x) * (d_x - point1_x) \
                + (point2_y - point1_y) * (d_y - point1_y)
        dist2 = (point2_x - point1_x) ** 2 + (point2_y - point1_y) ** 2

        if cross <= 0:
            return math.sqrt((d_x - point1_x) ** 2 + (d_y - point1_y) ** 2)
        if cross >= dist2:
            return math.sqrt((d_x - point2_x) ** 2 + (d_y - point2_y) ** 2)
        r = cross / dist2
        p_x = point1_x + (point2_x - point1_x) * r
        p_y = point1_y + (point2_y - point1_y) * r
        return math.sqrt((d_x - p_x) ** 2 + (d_y - p_y) ** 2)

    def reorganize(self, p, d):
        """
        Reorganize matrix P to canonical form by grouping
        absorbing states into ABS and transient states into TR.

        :param np.array p: matrix P
        :param Point d: the destination node
        :return: matrix Q(TR * TR), matrix S(TR * ABS)

        """
        ABS = []
        TR = []
        absorbing_state = [-1 for j in range(len(self.edges))]
        for index, node in enumerate(self.nodes):
            if node == d or self.edges[index] == absorbing_state:
                ABS.append(index)
            else:
                TR.append(index)

        p_left_top = p[np.ix_(TR, TR)]
        p_left_bottom = p[np.ix_(ABS, TR)]
        p_right_top = p[np.ix_(TR, ABS)]
        p_right_bottom = p[np.ix_(ABS, ABS)]

        return p_left_top, p_right_top

    def matrix_multiply(self, a, n):
        """
        Calculate n power of matrix A.

        :param np.array a: matrix A
        :param int n: n power
        :return: matrix result

        """
        result = np.identity(a.shape[0])
        for i in range(n):
            result = np.dot(result, a)
        return result

    def step_t(self):
        """
        Set step t as the diameter of the transfer network by Floyd.

        :return: the diameter

        """
        edge_matrix_len = len(self.edges)
        weight = [[sys.maxsize for j in range(edge_matrix_len)]
                  for i in range(edge_matrix_len)]
        for i in range(edge_matrix_len):
            weight[i][i] = 0
            for j in range(edge_matrix_len):
                if self.edges[i][j] != -1 and i != j:
                    weight[i][j] = 1

        for k in range(edge_matrix_len):
            for i in range(edge_matrix_len):
                for j in range(edge_matrix_len):
                    if weight[i][j] > weight[i][k] + weight[k][j]:
                        weight[i][j] = weight[i][k] + weight[k][j]
        return max(max(weight))

    def cal_vector(self, d, p, q):
        """
        Get the column vector V of each node through matrix P and matrix Q.

        :param Point d: the node d
        :param np.array p: matrix P
        :param np.array q: matrix Q
        :return: column vector V

        """
        TR = []
        absorbing_state = [-1 for j in range(len(self.edges))]

        # D=S[*,d]
        for index, node in enumerate(self.nodes):
            if not (node == d or self.edges[index] == absorbing_state):
                TR.append(index)
        D = p[np.ix_(TR, [self.nodes.index(d)])]

        # V=D+Q·D+Q^2·D+...+Q^(t-1)·D
        v = np.zeros(D.shape)
        for j in range(0, self.step_t()):
            v = v + np.dot(self.matrix_multiply(q, j), D)
        return v

max_stars_repo_name: Clloud/MostPopularRoute

max_issues_repo_name: Clloud/popular-route

max_forks_repo_name: Clloud/popular-route