Hot answers tagged

31

It sounds like you want to simulate from a truncated distribution, and in your specific example, a truncated normal. There are a variety of methods for doing so, some simple, some relatively efficient. I'll illustrate some approaches on your normal example. Here's one very simple method for generating one at a time (in some kind of pseudocode): $\tt{...


22

What you are doing is wrong: it does not make sense to compute PRESS for PCA like that! Specifically, the problem lies in your step #5. Naïve approach to PRESS for PCA Let the data set consist of $n$ points in $d$-dimensional space: $\mathbf x^{(i)} \in \mathbb R^d, \, i=1 \dots n$. To compute reconstruction error for a single test data point $\mathbf x^{(...


21

Nominal vs Interval The most classic "correlation" measure between a nominal and an interval ("numeric") variable is Eta, also called correlation ratio, and equal to the root R-square of the one-way ANOVA (with p-value = that of the ANOVA). Eta can be seen as a symmetric association measure, like correlation, because Eta of ANOVA (with the nominal as ...


20

The result that $p$ values have a uniform distribution under $H_0$ holds for continuously distributed test statistics - at least for point nulls, as you have here. As James Stanley mentions in comments the distribution of the test statistic is discrete, so that result doesn't apply. You may have no errors at all in your code (though I wouldn't display a ...


18

Should it be a 45 degree line? It depends! A QQ plot is the parametric curve defined by: \begin{align*} x &= F^{-1}(p)\\ y &= G^{-1}(p) \end{align*} for $p \in [0, 1]$. Where $F^{-1}$ and $G^{-1}$ are inverse CDF functions. If $F = G$ then $x(p)=y(p)$ and it would be on a 45 degree line. Another case... Let $\Phi^{-1}(p)$ be the standard normal ...


17

Assuming you have the original data and not just the summary of the fits, the general solution to this problem is to fit a model with an interaction, i.e. to go back to the data and fit the model $$ Y = \beta_0 + \beta_1 I(t>t_I) + \beta_2 (t-t_I) + \beta_3 I(t>t_I) (t-t_I) $$ where $I(t>t_I)$ is an indicator variable, i.e. =1 if $t>t_I$ and 0 ...


14

A kernel density estimator (KDE) produces a distribution that is a location mixture of the kernel distribution, so to draw a value from the kernel density estimate all you need do is (1) draw a value from the kernel density and then (2) independently select one of the data points at random and add its value to the result of (1). Here is the result of this ...


13

According to wikipedia, mutual information of two random variables may be calculated using the following formula: $$ I(X;Y) = \sum_{y \in Y} \sum_{x \in X} p(x,y) \log{ \left(\frac{p(x,y)}{p(x)\,p(y)} \right) } $$ If I pick up your code from this: [co1, ce1] = hist(randpoints1, bins); [co2, ce2] = hist(...


13

You're referring to a transformation from a pair of independent variates $(X,Y)$ to the polar representation $(R,\theta)$ (radius and angle), and then looking at the marginal distribution of $\theta$. I'm going to offer a somewhat intuitive explanation (though a mathematical derivation of the density does essentially what I describe informally). Note that ...


12

The facts that you are getting different answers from forward and backward selection, and that you get different answers when you change the seed, should give you pause. Clearly, these can't all be right. Most likely, none of them are. The simplest answer is that you should not use these methods at all. Here are some threads you might want to read: ...


11

Generating a PR curve is similar to generating an ROC curve. To draw such plots you need a full ranking of the test set. To make this ranking, you need a classifier which outputs a decision value rather than a binary answer. The decision value is a measure of confidence in a prediction which we can use to rank all test instances. As an example, the decision ...


10

Here is my implementation in R x <- c(1,2,1,1,3,4,4,1,2,4,1,4,3,4,4,4,3,1,3,2,3,3,3,4,2,2,3) xChar<-as.character(x) library(markovchain) mcX<-markovchainFit(xChar)$estimate mcX


10

ks.test in R allows one to adjust the mean and sd of the distribution to be tested against. e.g. x <- rnorm(1000, 4, 10) ks.test(x, "pnorm", mean = 4, sd = 10)


10

One prebuilt tool for visualizing high dimensional data is ggobi. It lets you color the points to represent groups and then has a few options for reducing the high dimensions to a 2 dimensional representation. One particularly nice tool is the 2D grand tour that basically rotates the data cloud in multiple dimensions and shows you an animation of the 2D ...


10

MATLAB was primarily developed for optimization and mathematical simulations in engineering problems. But yes it has performance issues when it comes to machine learning, optimization, etc., in terms of customizing ability. Over time, most statistical analysis / machine learning has shifted to R and Python, because of an active community presence for ...


10

Say your predictor matrix is $X$ and your response vector is $y$. PCA is concerned only with the (co)variance within the predictor matrix $X$ itself, while a regression model is (also) concerned with the covariance between $X$ and the response $y$. If there is no relationship between these concepts, dimension reduction by PCA can be harmful to your ...


9

You have interpreted these results correctly according to the conventional textbook scheme. Personally, I am often not a fan of the standard way of thinking about p-values. (Mounting soapbox...) Firstly, it's worth considering that there are several valid ways to look at p-values. Fisher thought of them as a continuous measure of evidence against the ...


9

Everything can be done with polyfit function only. Read about its optional output parameters in http://www.mathworks.nl/help/matlab/ref/polyfit.html For instance: [p,S,mu] = polyfit(x,y,n) where mu is the two-element vector [μ1,μ2], where μ1=mean(x), μ2=std(x) To compute error, you have to use another function taking output of polyfit: [y,delta] = ...


9

Remember that, if $Y_i$ are independent $\mathrm{Gamma}(a_i,b)$, for $i=1,\dots,k$, then $$ (X_1,\dots,X_k) = \left(\frac{Y_1}{\sum_{j=1}^k Y_j}, \dots, \frac{Y_k}{\sum_{j=1}^k Y_j} \right) \sim \mathrm{Dirichlet}(a_1,\dots,a_k) \, .$$ The proof can be found on page 594 of Luc Devroye's book. Therefore, one possibility is to compute a Monte Carlo ...


9

Here is the general (semi-parametric-bootstrap) algorithm in more detail: $\text{B}$ = number of bootstraps the model: $y = x\beta + \epsilon$ let $\hat{\epsilon}$ be the residuals Run the regression and obtain the estimator(s) $\hat\beta$ and residuals $\hat\epsilon$. Resample the residuals with replacement and obtain the bootstrapped residual vector $\...


9

This is based on the standard approximation to the Hessian of a nonlinear least squares problem used by Gauss-Newton and Levenberg-Marquardt algorithms. Consider the nonlinear least squares problem: minimize $1/2r(x)^Tr(x)$. Let $J$ = Jacobian of r(x). The Hessian of the objective = $J^TJ +$ higher order terms. The Gauss-Newton or Levenberg-Marquardt ...


9

Note that numpy.cov() considers its input data matrix to have observations in each column, and variables in each row, so to get numpy.cov() to return what other packages do, you have to pass the transpose of the data matrix to numpy.cov(). The Python code that you linked can be used to simulate what other packages do, but it contains some errors: N should ...


8

In R, you can just use the function ks.test with the following arguments: ks.test(your_data, "pnorm", mean=test_mu, sd=test_sd) Where your_data is your data vector, test_mu is the specific mean of the theoretical normal distribution and test_sd its standard deviation. To inspect your data graphically, you can use the function qqPlot from the car package. ...


8

You could give tSNE a try. It is pretty straightforward to use. It works with Octave, in addition to Matlab and Python. Take a look at the guide to get a first plot within a minute.


8

Here are a few thoughts: If you have sufficient data, you may not need your data to actually be normal, because the Central Limit Theorem will cover for you. (What counts as "sufficient data" will depend on the manner in which your data are non-normal, though.) You can test for a difference in the medians (rather than means) using the Mann-Whitney U-...


8

This will always happen for large-enough values of $\lambda$ (the regularization coefficient). (If your predictors aren't very good, this will be a very small value.) Are you providing it explicitly? The linked documentation says for the parameter Lambda: Vector of nonnegative Lambda values. See Definitions. – If you do not supply Lambda, lasso ...


8

This means, that optimization algorithm detected that with high probability (not in the strict, mathematical sense) you can speed up your training by turning the -h 0 flag in your options. Basically, -h is the shrinking heuristics, implemented in the libsvm package which for some data significantly reduces number of required computations, while in others - ...


8

You would use the Bonferroni for a one-way test. But let's be clear: You would not use the Bonferroni adjustment on the Kruskal-Wallis test itself. The Kruskal-Wallis test is an omnibus test, controlling for an overall false-positive rate. You would use the Bonferroni for post hoc Dunn's pairwise tests. Indeed, Dunn introduced the "Bonferroni" adjustment. ...


8

MATLAB always uses the centred and scaled variables for the computations within ridge. It just back-transforms them before returning them. As you have a really small matrix this probably makes a noticeable difference. You can reproduce the Python results in MATLAB easily: X = [1 1 2 ; 3 4 2 ; 6 5 2 ; 5 5 3]; Y = [1 0 0 1]; k = 10; % which is the ridge ...


8

That is correct because $b_{OLS}$ is the minimizer of MSE by definition. The problem ($X^TX$ is invertible here) has only one minimum and any value other than $b_{OLS}$ will have higher MSE on the training dataset.


Only top voted, non community-wiki answers of a minimum length are eligible