Compute fit residuals
184722
def fun_residuals(params, xnor, ynor, w, bbox, k, ext):
    """"""

    spl = LSQUnivariateSpline(
        x=xnor,
        y=ynor,
        t=[item.value for item in params.values()],
        w=w,
        bbox=bbox,
        k=k,
        ext=ext,
        check_finite=False
    )
    return spl.get_residual()
360446
def _bottleneck_residual(self, x, in_filter, out_filter, stride,
                           activate_before_residual=False):
    """"""
    if activate_before_residual:
      with tf.variable_scope('common_bn_relu'):
        x = self._layer_norm('init_bn', x)
        x = self._relu(x, self.hps.relu_leakiness)
        orig_x = x
    else:
      with tf.variable_scope('residual_bn_relu'):
        orig_x = x
        x = self._layer_norm('init_bn', x)
        x = self._relu(x, self.hps.relu_leakiness)

    with tf.variable_scope('sub1'):
      x = self._conv('conv1', x, 1, in_filter, out_filter / 4, stride)

    with tf.variable_scope('sub2'):
      x = self._layer_norm('bn2', x)
      x = self._relu(x, self.hps.relu_leakiness)
      x = self._conv('conv2', x, 3, out_filter / 4,
                     out_filter / 4, [1, 1, 1, 1])

    with tf.variable_scope('sub3'):
      x = self._layer_norm('bn3', x)
      x = self._relu(x, self.hps.relu_leakiness)
      x = self._conv('conv3', x, 1, out_filter /
                     4, out_filter, [1, 1, 1, 1])

    with tf.variable_scope('sub_add'):
      if in_filter != out_filter:
        orig_x = self._conv('project', orig_x, 1,
                            in_filter, out_filter, stride)
      x += orig_x

    return x
26084
def _fit_transform(self, X):
        
        X = self._matrix_adjust(X)

        X = check_array(
            X,
            accept_sparse='csc',
            force_all_finite=False,
            dtype=int
        )

        if X.min() < 0:
            raise ValueError("X needs to contain only non-negative integers.")

        n_samples, n_features = X.shape

        # Remember which values should not be replaced by the value 'other'
        if self.minimum_fraction is not None:
            do_not_replace_by_other = list()
            for column in range(X.shape[1]):
                do_not_replace_by_other.append(list())

                if sparse.issparse(X):
                    indptr_start = X.indptr[column]
                    indptr_end = X.indptr[column + 1]
                    unique = np.unique(X.data[indptr_start:indptr_end])
                    colsize = indptr_end - indptr_start
                else:
                    unique = np.unique(X[:, column])
                    colsize = X.shape[0]

                for unique_value in unique:
                    if np.isfinite(unique_value):
                        if sparse.issparse(X):
                            indptr_start = X.indptr[column]
                            indptr_end = X.indptr[column + 1]
                            count = np.nansum(unique_value ==
                                              X.data[indptr_start:indptr_end])
                        else:
                            count = np.nansum(unique_value == X[:, column])
                    else:
                        if sparse.issparse(X):
                            indptr_start = X.indptr[column]
                            indptr_end = X.indptr[column + 1]
                            count = np.nansum(~np.isfinite(
                                X.data[indptr_start:indptr_end]))
                        else:
                            count = np.nansum(~np.isfinite(X[:, column]))

                    fraction = float(count) / colsize
                    if fraction >= self.minimum_fraction:
                        do_not_replace_by_other[-1].append(unique_value)

                for unique_value in unique:
                    if unique_value not in do_not_replace_by_other[-1]:
                        if sparse.issparse(X):
                            indptr_start = X.indptr[column]
                            indptr_end = X.indptr[column + 1]
                            X.data[indptr_start:indptr_end][
                                X.data[indptr_start:indptr_end] ==
                                unique_value] = SPARSE_ENCODINGS['OTHER']
                        else:
                            X[:, column][X[:, column] == unique_value] = SPARSE_ENCODINGS['OTHER']

            self.do_not_replace_by_other_ = do_not_replace_by_other

        if sparse.issparse(X):
            n_values = X.max(axis=0).toarray().flatten() + len(SPARSE_ENCODINGS)
        else:
            n_values = np.max(X, axis=0) + len(SPARSE_ENCODINGS)

        self.n_values_ = n_values
        n_values = np.hstack([[0], n_values])
        indices = np.cumsum(n_values)
        self.feature_indices_ = indices

        if sparse.issparse(X):
            row_indices = X.indices
            column_indices = []
            for i in range(len(X.indptr) - 1):
                nbr = X.indptr[i+1] - X.indptr[i]
                column_indices_ = [indices[i]] * nbr
                column_indices_ += X.data[X.indptr[i]:X.indptr[i+1]]
                column_indices.extend(column_indices_)
            data = np.ones(X.data.size)
        else:
            column_indices = (X + indices[:-1]).ravel()
            row_indices = np.repeat(np.arange(n_samples, dtype=np.int32),
                                    n_features)
            data = np.ones(n_samples * n_features)

        out = sparse.coo_matrix((data, (row_indices, column_indices)),
                                shape=(n_samples, indices[-1]),
                                dtype=self.dtype).tocsc()

        mask = np.array(out.sum(axis=0)).ravel() != 0
        active_features = np.where(mask)[0]
        out = out[:, active_features]
        self.active_features_ = active_features
        return out.tocsr() if self.sparse else out.toarray()
542179
def fit(self, X, *args, **kwargs):
        

        self.constant_value = self._get_constant_value(X)

        if self.constant_value is None:
            if self.unfittable_model:
                self.model = getattr(scipy.stats, self.model_class)(*args, **kwargs)
            else:
                self.model = getattr(scipy.stats, self.model_class)(X, *args, **kwargs)

            for name in self.METHOD_NAMES:
                attribute = getattr(self.__class__, name)
                if isinstance(attribute, str):
                    setattr(self, name, getattr(self.model, attribute))

                elif attribute is None:
                    setattr(self, name, missing_method_scipy_wrapper(lambda x: x))

        else:
            self._replace_constant_methods()

        self.fitted = True
1139955
def fit(self):
        """"""
        self._mcmcfit = self.mcmcsetup.run()
        self._mcmcfit.burnin(self.burnin)
        dmin = min(self._mcmcfit.depth_segments)
        dmax = max(self._mcmcfit.depth_segments)
        self._thick = (dmax - dmin) / len(self.mcmcfit.depth_segments)
        self._depth = np.arange(dmin, dmax + 0.001)
        self._age_ensemble = np.array([self.agedepth(d=dx) for dx in self.depth])
383186
def fit_transform(self, X, y=None, **fit_params):
        
        y = column_or_1d(X, warn=True)
        _check_numpy_unicode_bug(X)
        self.classes_, X = np.unique(X, return_inverse=True)
        return X
47670
def _evaluate_fits(self, x, fit_key):
        
        fit = self.fits[fit_key]
        fit_params = self._get_fit_params(np.copy(x), fit_key)
        if type(fit) == list:
            res = []
            for i in range(len(fit)):
                res.append(fit[i](fit_params))
            return np.array(res)
        else:
            return fit(fit_params)
45113
def fitlin_clipped(xy,uv,verbose=False,mode='rscale',nclip=3,reject=3):
    
    fitting_funcs = {'rscale':fitlin_rscale,'general':fitlin}
    # Get the fitting function to be used
    fit_func = fitting_funcs[mode.lower()]
    # Perform the initial fit
    P,Q = fit_func(xy,uv)
    xyc = apply_fitlin(xy,P,Q)

    # compute residuals from fit for input positions
    dx = uv[:,0] - xyc[0]
    dy = uv[:,1] - xyc[1]
    fit_rms = [dx.std(),dy.std()]

    if nclip > 0:
        data = xy.copy()
        outdata = uv.copy()

    numclipped = 0
    for i in range(nclip):
        iterclipped = 0
        xyc = apply_fitlin(data,P,Q)

        # compute residuals from fit for input positions
        dx = outdata[:,0] - xyc[0]
        dy = outdata[:,1] - xyc[1]

        # find indices of outliers in x and y
        xout = np.where(np.abs(dx - dx.mean()) > reject*dx.std())
        yout = np.where(np.abs(dy - dy.mean()) > reject*dy.std())
        # concatenate those indices and sort them
        outliers_indx = xout[0].tolist()+yout[0].tolist()
        outliers_indx.sort()
        # define the full range of indices for the data points left
        full_indx = list(range(data.shape[0]))
        # remove all unique indices specified in outliers from full range
        for o in outliers_indx:
            # only remove if it has not been removed already
            # accounts for the same point being an outlier in both x and y
            if full_indx.count(o) > 0:
                full_indx.remove(o)
            iterclipped += 1

        if iterclipped == 0:
            break

        numclipped += iterclipped
        if verbose:
            print('Removed a total of ',numclipped,' points through iteration ',i+1)
        # create clipped data
        data_iter = np.zeros([len(full_indx),2],dtype=data.dtype)
        if verbose:
            print('Iter #',i+1,' data:',data.shape,data_iter.shape,len(full_indx))

        data_iter[:,0] = data[:,0][full_indx]
        data_iter[:,1] = data[:,1][full_indx]
        outdata_iter = np.zeros([len(full_indx),2],dtype=data.dtype)
        outdata_iter[:,0] = outdata[:,0][full_indx]
        outdata_iter[:,1] = outdata[:,1][full_indx]

        # perform the fit again with the clipped data and go to the next iteration
        data = data_iter
        outdata = outdata_iter
        P,Q = fit_func(data,outdata)
        # compute residuals from fit for input positions
        xyc = apply_fitlin(data,P,Q)
        dx = outdata[:,0] - xyc[0]
        dy = outdata[:,1] - xyc[1]
        fit_rms = [dx.std(),dy.std()]

    if verbose:
        print('Fit clipped ',numclipped,' points over ',nclip,' iterations.')
    return P,Q,fit_rms
246129
def residual_norm(A, x, b):
    """"""
    return norm(np.ravel(b) - A*np.ravel(x))
497489
def fit_transform(self,
                      X,
                      num_epochs=10,
                      updates_epoch=10,
                      stop_param_updates=dict(),
                      batch_size=1,
                      show_progressbar=False,
                      show_epoch=False):
        """"""
        self.fit(X,
                 num_epochs,
                 updates_epoch,
                 stop_param_updates,
                 batch_size,
                 show_progressbar,
                 show_epoch)

        return self.transform(X, batch_size=batch_size)
931908
public void fit(INDArray x) {
        this.data = x;
        for(int i = 0; i < numTrees; i++) {
            RPTree tree = new RPTree(data.columns(),maxSize,similarityFunction);
            tree.buildTree(x);
            trees.add(tree);
        }
    }
879837
def _fit_stage_componentwise(X, residuals, sample_weight, **fit_params):
    """"""
    n_features = X.shape[1]

    base_learners = []
    error = numpy.empty(n_features)
    for component in range(n_features):
        learner = ComponentwiseLeastSquares(component).fit(X, residuals, sample_weight)
        l_pred = learner.predict(X)
        error[component] = squared_norm(residuals - l_pred)
        base_learners.append(learner)

    # TODO: could use bottleneck.nanargmin for speed
    best_component = numpy.nanargmin(error)
    best_learner = base_learners[best_component]
    return best_learner
758851
def _compute_delta_beta(self, X, T, E, weights, index=None):
        
        score_residuals = self._compute_score(X, T, E, weights, index=index)

        d = X.shape[1]
        scaled_variance_matrix = self.variance_matrix_ * np.tile(self._norm_std.values, (d, 1)).T

        delta_betas = score_residuals.dot(scaled_variance_matrix)
        delta_betas.columns = self.hazards_.index

        return delta_betas
797140
def compute_resids(xy,uv,fit):
    
    print('FIT coeffs: ',fit['coeffs'])
    xn,yn = apply_fit(uv,fit['coeffs'])
    resids = xy - np.transpose([xn,yn])
    return resids
413876
def plot_series_residuals_internal(self, varied_data, varied_idx, **kwargs):
        
        nf = len(self.f_cb(*self.pre_process(
            self.internal_xout[0], self.internal_params_out[0])))
        xerr = np.empty((self.internal_xout.shape[0], nf))
        for idx, (res, params) in enumerate(zip(self.internal_xout, self.internal_params_out)):
            xerr[idx, :] = self.f_cb(res, params)
        return self.plot_series(xerr, varied_data, varied_idx, **kwargs)
54615
def _fit_diagV_noGP(
            self, XTY, XTDY, XTFY, YTY_diag, YTDY_diag, YTFY_diag,
            XTX, XTDX, XTFX, X, Y, X_base, X_res, D, F, run_TRs,
            current_vec_U_chlsk_l,
            current_a1, current_logSNR2,
            idx_param_fitU, idx_param_fitV,
            l_idx, n_C, n_T, n_V, n_l, n_run, n_X0, rank):
        
        init_iter = self.init_iter
        logger.info('second fitting without GP prior'
                    ' for {} times'.format(init_iter))

        # Initial parameters
        param0_fitU = np.empty(
            np.sum(np.size(v) for v in idx_param_fitU.values()))
        param0_fitV = np.empty(np.size(idx_param_fitV['log_SNR2']))
        # We cannot use the same logic as the line above because
        # idx_param_fitV also includes entries for GP parameters.
        param0_fitU[idx_param_fitU['Cholesky']] = \
            current_vec_U_chlsk_l.copy()
        param0_fitU[idx_param_fitU['a1']] = current_a1.copy()
        param0_fitV[idx_param_fitV['log_SNR2']] = \
            current_logSNR2[:-1].copy()

        L = np.zeros((n_C, rank))
        tol = self.tol * 5
        for it in range(0, init_iter):
            X0TX0, X0TDX0, X0TFX0, XTX0, XTDX0, XTFX0, \
                X0TY, X0TDY, X0TFY, X0, X_base, n_X0, _ \
                = self._prepare_data_XYX0(
                    X, Y, X_base, X_res, D, F, run_TRs, no_DC=True)

            # fit U, the covariance matrix, together with AR(1) param
            param0_fitU[idx_param_fitU['Cholesky']] = \
                current_vec_U_chlsk_l \
                + self.random_state_.randn(n_l) \
                * np.linalg.norm(current_vec_U_chlsk_l) \
                / n_l**0.5 * np.exp(-it / init_iter * self.anneal_speed - 1)
            param0_fitU[idx_param_fitU['a1']] = current_a1
            res_fitU = scipy.optimize.minimize(
                self._loglike_AR1_diagV_fitU, param0_fitU,
                args=(XTX, XTDX, XTFX, YTY_diag, YTDY_diag, YTFY_diag,
                      XTY, XTDY, XTFY, X0TX0, X0TDX0, X0TFX0,
                      XTX0, XTDX0, XTFX0, X0TY, X0TDY, X0TFY,
                      current_logSNR2, l_idx, n_C,
                      n_T, n_V, n_run, n_X0, idx_param_fitU, rank),
                method=self.optimizer, jac=True, tol=tol,
                options=self.minimize_options)
            current_vec_U_chlsk_l = \
                res_fitU.x[idx_param_fitU['Cholesky']]
            current_a1 = res_fitU.x[idx_param_fitU['a1']]
            norm_fitUchange = np.linalg.norm(res_fitU.x - param0_fitU)
            logger.debug('norm of parameter change after fitting U: '
                         '{}'.format(norm_fitUchange))
            param0_fitU = res_fitU.x.copy()

            # fit V, reflected in the log(SNR^2) of each voxel
            rho1 = np.arctan(current_a1) * 2 / np.pi
            L[l_idx] = current_vec_U_chlsk_l
            X0TAX0, XTAX0, X0TAY, X0TAX0_i, \
                XTAcorrX, XTAcorrY, YTAcorrY, \
                LTXTAcorrY, XTAcorrXL, LTXTAcorrXL = \
                self._precompute_ar1_quad_forms(XTY, XTDY, XTFY,
                                                YTY_diag, YTDY_diag, YTFY_diag,
                                                XTX, XTDX, XTFX,
                                                X0TX0, X0TDX0, X0TFX0,
                                                XTX0, XTDX0, XTFX0,
                                                X0TY, X0TDY, X0TFY,
                                                L, rho1, n_V, n_X0)
            res_fitV = scipy.optimize.minimize(
                self._loglike_AR1_diagV_fitV, param0_fitV,
                args=(X0TAX0, XTAX0, X0TAY,
                      X0TAX0_i, XTAcorrX, XTAcorrY, YTAcorrY,
                      LTXTAcorrY, XTAcorrXL, LTXTAcorrXL,
                      current_vec_U_chlsk_l,
                      current_a1, l_idx, n_C, n_T, n_V, n_run,
                      n_X0, idx_param_fitV, rank,
                      False, False),
                method=self.optimizer, jac=True, tol=tol,
                options=self.minimize_options)

            current_logSNR2[0:n_V - 1] = res_fitV.x
            current_logSNR2[-1] = - np.sum(current_logSNR2[0:n_V - 1])

            norm_fitVchange = np.linalg.norm(res_fitV.x - param0_fitV)
            logger.debug('norm of parameter change after fitting V: '
                         '{}'.format(norm_fitVchange))
            logger.debug('E[log(SNR2)^2]: {}'.format(
                np.mean(current_logSNR2**2)))

            # The lines below are for debugging purpose.
            # If any voxel's log(SNR^2) gets to non-finite number,
            # something might be wrong -- could be that the data has
            # nothing to do with the design matrix.
            if np.any(np.logical_not(np.isfinite(current_logSNR2))):
                logger.warning('Initial fitting: iteration {}'.format(it))
                logger.warning('current log(SNR^2): '
                               '{}'.format(current_logSNR2))
                logger.warning('log(sigma^2) has non-finite number')

            param0_fitV = res_fitV.x.copy()

            # Re-estimating X_res from residuals
            current_SNR2 = np.exp(current_logSNR2)
            if self.auto_nuisance:
                LL, LAMBDA_i, LAMBDA, YTAcorrXL_LAMBDA, current_sigma2 \
                    = self._calc_LL(rho1, LTXTAcorrXL, LTXTAcorrY, YTAcorrY,
                                    X0TAX0, current_SNR2,
                                    n_V, n_T, n_run, rank, n_X0)
                betas = current_SNR2 * np.dot(L, YTAcorrXL_LAMBDA.T)
                beta0s = np.einsum(
                    'ijk,ki->ji', X0TAX0_i,
                    (X0TAY - np.einsum('ikj,ki->ji', XTAX0, betas)))
                residuals = Y - np.dot(X, betas) - np.dot(
                    X_base, beta0s[:np.shape(X_base)[1], :])
                X_res = self.nureg_method(
                    self.n_nureg_).fit_transform(
                    self.preprocess_residual(residuals))
            if norm_fitVchange / np.sqrt(param0_fitV.size) < tol \
                    and norm_fitUchange / np.sqrt(param0_fitU.size) \
                    < tol:
                break
        return current_vec_U_chlsk_l, current_a1, current_logSNR2, X_res
41025
def _compute_errors(self):
        

        try:
            coeffs = fit_first_and_second_harmonics(self.sample.values[0],
                                                    self.sample.values[2])
            covariance = coeffs[1]
            coeffs = coeffs[0]
            model = first_and_second_harmonic_function(self.sample.values[0],
                                                       coeffs)
            residual_rms = np.std(self.sample.values[2] - model)
            errors = np.diagonal(covariance) * residual_rms

            eps = self.sample.geometry.eps
            pa = self.sample.geometry.pa

            # parameter errors result from direct projection of
            # coefficient errors.  These showed to be the error estimators
            # that best convey the errors measured in Monte Carlo
            # experiments (see Busko 1996; ASPC 101, 139).
            ea = abs(errors[2] / self.grad)
            eb = abs(errors[1] * (1. - eps) / self.grad)
            self.x0_err = np.sqrt((ea * np.cos(pa))**2 + (eb * np.sin(pa))**2)
            self.y0_err = np.sqrt((ea * np.sin(pa))**2 + (eb * np.cos(pa))**2)
            self.ellip_err = (abs(2. * errors[4] * (1. - eps) / self.sma /
                                  self.grad))
            if (abs(eps) > np.finfo(float).resolution):
                self.pa_err = (abs(2. * errors[3] * (1. - eps) / self.sma /
                                   self.grad / (1. - (1. - eps)**2)))
            else:
                self.pa_err = 0.
        except Exception:    # we want to catch everything
            self.x0_err = self.y0_err = self.pa_err = self.ellip_err = 0.