Hard thresholding with projection onto non-negative numbers
976744
def prox_hard_plus(X, step, thresh=0):
    
    return prox_plus(prox_hard(X, step, thresh=thresh), step)
600779
def mask_unphysical(self, data):
        """"""
        if not self.valid_range:
            return data
        else:
            return np.ma.masked_outside(data, np.min(self.valid_range),
                                        np.max(self.valid_range))
1174634
private ApproximationLine approximateLowerHull(ConvexHull convexHull, double[] log_k, double sum_log_k, double sum_log_k2, double[] log_kDist, double sum_log_kDist, double sum_log_k_kDist) {
    // StringBuilder msg = new StringBuilder(1000);
    int[] lowerHull = convexHull.getLowerHull();
    int l = convexHull.getNumberOfPointsInLowerHull();
    int k_0 = settings.kmax - lowerHull.length + 1;

    // linear search on all line segments on the lower convex hull
    // msg.append("lower hull l = ").append(l).append('\n');
    double low_error = Double.MAX_VALUE;
    double low_m = 0.0;
    double low_t = 0.0;

    for(int i = 1; i < l; i++) {
      double cur_m = (log_kDist[lowerHull[i]] - log_kDist[lowerHull[i - 1]]) / (log_k[lowerHull[i]] - log_k[lowerHull[i - 1]]);
      double cur_t = log_kDist[lowerHull[i]] - cur_m * log_k[lowerHull[i]];
      double cur_error = ssqerr(k_0, settings.kmax, log_k, log_kDist, cur_m, cur_t);
      // msg.append(" Segment = ").append(i).append(" m =
      // ").append(cur_m).append(" t = ").append(cur_t).append(" lowerror =
      // ").append(cur_error).append('\n');
      if(cur_error < low_error) {
        low_error = cur_error;
        low_m = cur_m;
        low_t = cur_t;
      }
    }

    // linear search on all points of the lower convex hull
    boolean is_right = true; // NEEDED FOR PROOF CHECK
    for(int i = 0; i < l; i++) {
      double cur_m = optimize(k_0, settings.kmax, sum_log_k, sum_log_k2, log_k[lowerHull[i]], log_kDist[lowerHull[i]], sum_log_k_kDist, sum_log_kDist);
      double cur_t = log_kDist[lowerHull[i]] - cur_m * log_k[lowerHull[i]];
      // only valid if both neighboring points are underneath y=mx+t
      if((i == 0 || log_kDist[lowerHull[i - 1]] >= log_kDist[lowerHull[i]] - cur_m * (log_k[lowerHull[i]] - log_k[lowerHull[i - 1]])) && (i == l - 1 || log_kDist[lowerHull[i + 1]] >= log_kDist[lowerHull[i]] + cur_m * (log_k[lowerHull[i + 1]] - log_k[lowerHull[i]]))) {
        double cur_error = ssqerr(k_0, settings.kmax, log_k, log_kDist, cur_m, cur_t);
        if(cur_error < low_error) {
          low_error = cur_error;
          low_m = cur_m;
          low_t = cur_t;
        }
      }

      // check proof of bisection search
      if(!(i > 0 && log_kDist[lowerHull[i - 1]] < log_kDist[lowerHull[i]] - cur_m * (log_k[lowerHull[i]] - log_k[lowerHull[i - 1]])) && !is_right) {
        LOG.warning("ERROR lower: The bisection search will not work properly!");
        if(!(i < l - 1 && log_kDist[lowerHull[i + 1]] < log_kDist[lowerHull[i]] + cur_m * (log_k[lowerHull[i + 1]] - log_k[lowerHull[i]]))) {
          is_right = false;
        }
      }
    }
    return new ApproximationLine(k_0, low_m, low_t);
  }
1042886
private static void thresholdBlock(byte[] luminances,
                                     int xoffset,
                                     int yoffset,
                                     int threshold,
                                     int stride,
                                     BitMatrix matrix) {
    for (int y = 0, offset = yoffset * stride + xoffset; y < BLOCK_SIZE; y++, offset += stride) {
      for (int x = 0; x < BLOCK_SIZE; x++) {
        // Comparison needs to be <= so that black == 0 pixels are black even if the threshold is 0.
        if ((luminances[offset + x] & 0xFF) <= threshold) {
          matrix.set(xoffset + x, yoffset + y);
        }
      }
    }
  }
225582
function decision2DStumpTrain(data, labels, ix, options) {
    
    options = options || {};
    var numtries = options.numTries || 10;
    
    // choose a dimension at random and pick a best split
    var N = ix.length;
    
    var ri1= 0;
    var ri2= 1;
    if(data[0].length > 2) {
      // more than 2D data. Pick 2 random dimensions
      ri1= randi(0, data[0].length);
      ri2= randi(0, data[0].length);
      while(ri2 == ri1) ri2= randi(0, data[0].length); // must be distinct!
    }
    
    // evaluate class entropy of incoming data
    var H= entropy(labels, ix);
    var bestGain=0; 
    var bestw1, bestw2, bestthr;
    var dots= new Array(ix.length);
    for(var i=0;i<numtries;i++) {
        
        // pick random line parameters
        var alpha= randf(0, 2*Math.PI);
        var w1= Math.cos(alpha);
        var w2= Math.sin(alpha);
        
        // project data on this line and get the dot products
        for(var j=0;j<ix.length;j++) {
          dots[j]= w1*data[ix[j]][ri1] + w2*data[ix[j]][ri2];
        }
        
        // we are in a tricky situation because data dot product distribution
        // can be skewed. So we don't want to select just randomly between
        // min and max. But we also don't want to sort as that is too expensive
        // let's pick two random points and make the threshold be somewhere between them.
        // for skewed datasets, the selected points will with relatively high likelihood
        // be in the high-desnity regions, so the thresholds will make sense
        var ix1= ix[randi(0, N)];
        var ix2= ix[randi(0, N)];
        while(ix2==ix1) ix2= ix[randi(0, N)]; // enforce distinctness of ix2
        var a= Math.random();
        var dotthr= dots[ix1]*a + dots[ix2]*(1-a);
        
        // measure information gain we'd get from split with thr
        var l1=1, r1=1, lm1=1, rm1=1; //counts for Left and label 1, right and label 1, left and minus 1, right and minus 1
        for(var j=0;j<ix.length;j++) {
            if(dots[j] < dotthr) {
              if(labels[ix[j]]==1) l1++;
              else lm1++;
            } else {
              if(labels[ix[j]]==1) r1++;
              else rm1++;
            }
        }
        var t= l1+lm1; 
        l1=l1/t;
        lm1=lm1/t;
        t= r1+rm1;
        r1=r1/t;
        rm1= rm1/t;
        
        var LH= -l1*Math.log(l1) -lm1*Math.log(lm1); // left and right entropy
        var RH= -r1*Math.log(r1) -rm1*Math.log(rm1);
        
        var informationGain= H - LH - RH;
        //console.log("Considering split %f, entropy %f -> %f, %f. Gain %f", thr, H, LH, RH, informationGain);
        if(informationGain > bestGain || i === 0) {
            bestGain= informationGain;
            bestw1= w1;
            bestw2= w2;
            bestthr= dotthr;
        }
    }
    
    model= {};
    model.w1= bestw1;
    model.w2= bestw2;
    model.dotthr= bestthr;
    return model;
}
1335060
def single_discriminator(x, filters=128, kernel_size=8,
                         strides=4, pure_mean=False):
  """"""
  with tf.variable_scope("discriminator"):
    net = layers().Conv2D(
        filters, kernel_size, strides=strides, padding="SAME", name="conv1")(x)
    if pure_mean:
      net = tf.reduce_mean(net, [1, 2])
    else:
      net = mean_with_attention(net, "mean_with_attention")
    return net
221119
def neighborhood_probability(self, threshold, radius):
        
        weights = disk(radius, dtype=np.uint8)
        thresh_data = np.zeros(self.data.shape[1:], dtype=np.uint8)
        neighbor_prob = np.zeros(self.data.shape, dtype=np.float32)
        for t in np.arange(self.data.shape[0]):
            thresh_data[self.data[t] >= threshold] = 1
            maximized = fftconvolve(thresh_data, weights, mode="same")
            maximized[maximized > 1] = 1
            maximized[maximized < 1] = 0
            neighbor_prob[t] = fftconvolve(maximized, weights, mode="same")
            thresh_data[:] = 0
        neighbor_prob[neighbor_prob < 1] = 0
        neighbor_prob /= weights.sum()
        return neighbor_prob
1144573
def pct_negative(self, threshold=0.0):
        

        return np.count_nonzero(self[self < threshold]) / self.count()
1160069
def estimateFromImages(imgs1, imgs2=None, mn_mx=None, nbins=100):
    '''
    
    '''
    if imgs2 is None:
        imgs2 = [None] * len(imgs1)
    else:
        assert len(imgs1) == len(imgs2)

    y_vals = np.empty((len(imgs1), nbins))
    w_vals = np.zeros((len(imgs1), nbins))

    if mn_mx is None:
        print('estimating min and max image value')
        mn = 1e6
        mx = -1e6
        # get min and max image value checking all first images:
        for n, i1 in enumerate(imgs1):
            print('%s/%s' % (n + 1, len(imgs1)))
            i1 = imread(i1)
            mmn, mmx = _getMinMax(i1)
            mn = min(mn, mmn)
            mx = mx = max(mx, mmx)
        print('--> min(%s), max(%s)' % (mn, mx))
    else:
        mn, mx = mn_mx

    x = None
    print('get noise level function')
    for n, (i1, i2) in enumerate(zip(imgs1, imgs2)):
        print('%s/%s' % (n + 1, len(imgs1)))

        i1 = imread(i1)
        if i2 is not None:
            i2 = imread(i2)

        x, y, weights, _ = calcNLF(i1, i2, mn_mx_nbins=(mn, mx, nbins), x=x)
        y_vals[n] = y
        w_vals[n] = weights

    # filter empty places:
    filledPos = np.sum(w_vals, axis=0) != 0
    w_vals = w_vals[:, filledPos]
    y_vals = y_vals[:, filledPos]
    x = x[filledPos]

    y_avg = np.average(np.nan_to_num(y_vals),
                       weights=w_vals,
                       axis=0)
    w_vals = np.sum(w_vals, axis=0)
    w_vals /= w_vals.sum()

    fitParams, fn, i = _evaluate(x, y_avg, w_vals)
    return x, fn, y_avg, y_vals, w_vals, fitParams, i
794716
function mask (contrast) {
	var a = Math.pow(392.498 * contrast, 0.7);
	var b = Math.pow(0.0153 * a, 4.0);
	return Math.pow(1.0 + b, 0.25);
}
38910
def generate_threshold_mask(hist):
    '''
    '''
    masked_array = np.ma.masked_values(hist, 0)
    masked_array = np.ma.masked_greater(masked_array, 10 * np.ma.median(hist))
    logging.info('Masking %d pixel(s)', np.ma.count_masked(masked_array))
    return np.ma.getmaskarray(masked_array)
1338767
def approximate_gradient(decision_function, sample, num_evals,
                         delta, constraint, shape, clip_min, clip_max):
  
  # Generate random vectors.
  noise_shape = [num_evals] + list(shape)
  if constraint == 'l2':
    rv = np.random.randn(*noise_shape)
  elif constraint == 'linf':
    rv = np.random.uniform(low=-1, high=1, size=noise_shape)

  axis = tuple(range(1, 1 + len(shape)))
  rv = rv / np.sqrt(np.sum(rv ** 2, axis=axis, keepdims=True))
  perturbed = sample + delta * rv
  perturbed = clip_image(perturbed, clip_min, clip_max)
  rv = (perturbed - sample) / delta

  # query the model.
  decisions = decision_function(perturbed)
  decision_shape = [len(decisions)] + [1] * len(shape)
  fval = 2 * decisions.astype(np_dtype).reshape(decision_shape) - 1.0

  # Baseline subtraction (when fval differs)
  if np.mean(fval) == 1.0:  # label changes.
    gradf = np.mean(rv, axis=0)
  elif np.mean(fval) == -1.0:  # label not change.
    gradf = - np.mean(rv, axis=0)
  else:
    fval = fval - np.mean(fval)
    gradf = np.mean(fval * rv, axis=0)

  # Get the gradient direction.
  gradf = gradf / np.linalg.norm(gradf)

  return gradf
1019621
def cuts_from_bbox(mask_nii, cuts=3):
    """"""
    from nibabel.affines import apply_affine

    mask_data = mask_nii.get_data() > 0.0

    # First, project the number of masked voxels on each axes
    ijk_counts = [
        mask_data.sum(2).sum(1),  # project sagittal planes to transverse (i) axis
        mask_data.sum(2).sum(0),  # project coronal planes to to longitudinal (j) axis
        mask_data.sum(1).sum(0),  # project axial planes to vertical (k) axis
    ]

    # If all voxels are masked in a slice (say that happens at k=10),
    # then the value for ijk_counts for the projection to k (ie. ijk_counts[2])
    # at that element of the orthogonal axes (ijk_counts[2][10]) is
    # the total number of voxels in that slice (ie. Ni x Nj).
    # Here we define some thresholds to consider the plane as "masked"
    # The thresholds vary because of the shape of the brain
    # I have manually found that for the axial view requiring 30%
    # of the slice elements to be masked drops almost empty boxes
    # in the mosaic of axial planes (and also addresses #281)
    ijk_th = [
        int((mask_data.shape[1] * mask_data.shape[2]) * 0.2),   # sagittal
        int((mask_data.shape[0] * mask_data.shape[2]) * 0.0),   # coronal
        int((mask_data.shape[0] * mask_data.shape[1]) * 0.3),   # axial
    ]

    vox_coords = []
    for ax, (c, th) in enumerate(zip(ijk_counts, ijk_th)):
        B = np.argwhere(c > th)
        if B.size:
            smin, smax = B.min(), B.max()

        # Avoid too narrow selections of cuts (very small masks)
        if not B.size or (th > 0 and (smin + cuts + 1) >= smax):
            B = np.argwhere(c > 0)

        # Resort to full plane if mask is seemingly empty
        smin, smax = B.min(), B.max() if B.size else (0, mask_data.shape[ax])
        inc = (smax - smin) / (cuts + 1)
        vox_coords.append([smin + (i + 1) * inc for i in range(cuts)])

    ras_coords = []
    for cross in np.array(vox_coords).T:
        ras_coords.append(apply_affine(
            mask_nii.affine, cross).tolist())
    ras_cuts = [list(coords) for coords in np.transpose(ras_coords)]
    return {k: v for k, v in zip(['x', 'y', 'z'], ras_cuts)}
883695
def _extract_peaks(specgram, neighborhood, threshold):
    
    kernel = np.ones(shape=neighborhood)
    local_averages = convolve(specgram, kernel / kernel.sum(), mode="constant", cval=0)

    # suppress all points below the floor value
    floor = (1 + threshold) * local_averages
    candidates = np.where(specgram > floor, specgram, 0)

    # grayscale dilation is equivalent to non-maximal suppression
    local_maximums = grey_dilation(candidates, footprint=kernel)
    peak_coords = np.argwhere(specgram == local_maximums)
    peaks = zip(peak_coords[:, 0], peak_coords[:, 1])

    return peaks
1236903
def norm(image):
    """"""
    min_value = image.min()
    max_value = image.max()
    if min_value == max_value:
        return image - min_value
    return (image - min_value) / (max_value - min_value)
1332315
private static ArrayList<double[]> getVisibleEdges(Projection2D proj, NumberVector s_min, NumberVector s_max) {
    final int dim = s_min.getDimensionality();
    double[] s_deltas = new double[dim];
    for(int i = 0; i < dim; i++) {
      s_deltas[i] = s_max.doubleValue(i) - s_min.doubleValue(i);
    }
    ArrayList<double[]> r_edges = new ArrayList<>(dim);
    for(int i = 0; i < dim; i++) {
      double[] delta = new double[dim];
      delta[i] = s_deltas[i];
      double[] deltas = proj.fastProjectRelativeDataToRenderSpace(delta);
      if(deltas[0] != 0 || deltas[1] != 0) {
        r_edges.add(deltas);
      }
    }
    return r_edges;
  }
1046981
def medianThreshold(img, threshold=0.1, size=3, condition='>', copy=True):
    '''
    
    '''
    from scipy.ndimage import median_filter

    indices = None
    if threshold > 0:
        blur = np.asfarray(median_filter(img, size=size))
        with np.errstate(divide='ignore', invalid='ignore', over='ignore'):

            if condition == '>':
                indices = abs((img - blur) / blur) > threshold
            else:
                indices = abs((img - blur) / blur) < threshold

        if copy:
            img = img.copy()

        img[indices] = blur[indices]
    return img, indices