Tag Info

$\mathbf{A}$ is an $(n+1) \times (n+1)$ matrix. It can be obtained as follows: $\textbf{A} = \left[ \matrix{1 & 1 & 1 & \cdots & 1 \cr x_0 & x_1 & x_2 & \cdots & x_{n} \cr \vdots & \vdots & \vdots & \cdots & \vdots \cr x_0^n & x_1^n & x_2^n & \cdots & ... 6 This doesn't randomly sample points, but instead chooses representative points deterministically. scipy.stats.norm.ppf(np.linspace(0, 1, 1000+2)[1:-1]) 6 It's a question of how you choose$\lambda$and$\Gamma$(of course). Think for a moment about what happens if you choose$\Gamma=I$and make$\lambda$large: in that case you say that it is more important to you to minimize the regularization term$\|Ix\|^2=\|x\|^2$than to minimize the misfit$\|Ax-b\|^2$. Obviously, making the term$\|x\|^2$small means ... 5 It's actually not all that complicated to calculate the Tracy-Widom CDF just from its definition: see On The Numerical Evaluation Of Fredholm Determinants by Folkmar Bornemann. The Wikipedia page gives the definition as $$F_2(s) = \det(I - A_s),$$ where$A_s$is the integral operator on$[s,\infty)$with the kernel $$\DeclareMathOperator{\Ai}{Ai} K(x,y) =... 5 Would a decomposition of the form A = XX^T suffice? This would be enough, e.g., if the end goal is sampling from the Gaussian distribution with this given covariance. If so, you can use the following formula, which is quite similar to your approximation:$$X = D^{1/2} + \frac{\sqrt{u^T D^{-1} u+1}-1}{u^T D^{-1} u} u u^T D^{-1/2}$$This follows from ... 4 Let \mathbf{\theta} be a Gaussian random vector with mean \mathbf{\mu} and covariance matrix \mathbf{\Sigma}_\mathbf{\theta}. Let \mathbf{p}_\theta denote the joint PDF. Let J_\mathbf{\theta} be the objective function, as its negative logarithm:$$J_\mathbf{\theta}=-ln(\mathbf{p}_\theta)$$By taking the partial derivatives w.r.t. \theta_d and ... 4 You can easily do this with for loops. Just use n to define your loop limits. Here's a fully functional solution for MATLAB (just define n and x first): %pre-allocate A A = zeros(n+1); %first row: for j=0:n A(1,j+1)=sum(x.^j); end %rows 2 through n for i=1:n A(i+1,1:n)=A(i,2:n+1); %copy from previous row A(i+1,n+1)=sum(x.^(n+i)); %compute last ... 4 Let f(\omega) be your power spectrum. Then maybe something like$$ \frac{\|f\|_{L^\infty}\|f\|_{L^0}}{\|f\|_{L^1}} = \frac{\mathrm{max}_{\omega\in\Omega} f(\omega)\cdot|\omega_{max}-\omega_{min}|}{\int_\Omega|f(\omega)|d\omega}. $$I know that L^0 isn't really good notation but I think it is useful for presenting this. This quantity is minimized by ... 3 You may use Kahan summation algorithm  The idea is to reschedule the sum operations in such a way precision loss is limited. The code is very simple (reproduced from  below). If this does not suffice, you may use multiprecision representations, such as quad doubles . They are supported by several languages / compilers (including GNU c). Finally, if ... 3 Generally, people refer to software packages as 'solver' that, well, solve some class of equations. For example, it could be stated that R contains a linear system solver. 3 Testing whether or not the mean is correct, or even if the histogram of your generated random variants "looks" like a certain distribution is not sufficient. Stick with much more rigorous test suites such as TestU01 or Diehard. Also, you really only have TWO random numbers in each row, because of the constraint that they sum to 1. This requires more ... 3 You don't need to sort your data first, if you have access to additional temporary storage, in which case, you should use a selection algorithm. 3 Mathematica has the TW distributions: http://reference.wolfram.com/language/ref/TracyWidomDistribution.html 3 In general, your type of question would be called a "multivariate goodness of fit test". If F(x_1,\ldots,x_n) is the n-dimensional CDF for the theoretical distribution, and the random variables (X_1,\ldots,X_n) are a sample from your ODE, then$$ Z_1 = F(X_1), \quad Z_2 = F(X_2\mid X_1), \quad\cdots\quad Z_n = F(X_n\mid X_1,\ldots,X_{n-1}) $$are ... 3 (Converted to an answer from my comments and expanded.) Basically you need to compute the fourth moments E[x_i x_j y_k y_l] for all i,j,k,l, given the second moments. These fourth moments are not universally determined, but they depend on the distribution of your variables, even in the one-dimensional case (that's exactly the reason why kurtosis is a ... 2 The AIC function need an 'lm' or 'glm' object (linear models). See functions lm and glm So just do : AIC(lm(logCPK~dataPOW)) 2 NumPy comes with a nifty random library with various distributions, including normal (Gaussian). From the Numpy documentation: mu, sigma = 0, 0.1 # mean and standard deviation s = np.random.normal(mu, sigma, 1000) which will give you 1000 normally distributed values with mean mu and standard deviation sigma. 2 Principal Components Analysis (PCA) is conducted using a Singular Value Decomposition (SVD) algorithm. As Bill Barth says above, the choice of sign of the principal component vectors is entirely arbitrary. 2 Established methodologies for benchmarking optimization software can be found in publications such as Benchmarking Optimization Software with Performance Profiles, Benchmarking Derivative-Free Optimization Algorithms, and Derivative-free optimization: a review of algorithms and comparison of software implementations. Generally speaking, algorithms are ... 2 You are adding positive terms, so you don't need to worry about the imprecision of the final result as if you had negative terms: (10^100+5)-(10^100+3)=2 but (10^100+5)+(10^100+3)=2*10^100 (as long you don't need 100 digits of precision). In your case I'd sort the exponents in decreasing order (in the sense of magnitude). Then I'd keep adding just in case ... 2 I think that what you are trying to do is to find a line passing through a set of data which is able to best fit that set of data. A Least Square approach, for instance, could be used for this purpose. Force a line to pass through the origin means that if you have the general affine equation for the line as$$ y = mx +q $$than the affine term q is null, ... 2 It's not entirely clear to me what step you struggle with. But for the sake of explanation, assume you have data x_i, i=1...N and you want to compute the normalized data y_i, i=1...N from the x_i. Then$$ y_i = \frac{x_i}{\sqrt{\sum_{j=1}^N x_j^2}}. $$The denominator is of course the same for every i, so you only have to compute it once. 2 Bagging, Boosting, and Bayesian Model Averaging/Combination are all widely used techniques for doing this. These are discussed in many textbooks on machine learning. 2 This is a typical use case for a paired t-test. The idea is to consider only the runtime difference \Delta t for each problem and test for the null hypothesis E(\Delta t)=0. For a step-by-step explanation, see e.g. (the article refers to segmentation evaluation, but on an abstract level the problem is identiclal to yours): Mao, Kanungo: "Empirical ... 2 What about the following? Dinv = np.diag(1.0 / np.sqrt(np.diag(cov_matrix))) corr = Dinv @ cov_matrix @ Dinv The above avoids any division except by the diagonal values, which should be nonzero anyway. 1 Mathematically, the probability density function (\operatorname{PDF}) for Z is given by the integral:$$\operatorname{PDF}(Z) = \int \delta\left(Z - f_a(X,Y)\right) \operatorname{PDF}(X,Y)\operatorname{d}X \operatorname{d}Y.$$If the transformation$$\left[\begin{array}{c} X \\ Y \end{array}\right] \rightarrow \left[\begin{array}{c} f_a(X,Y) \\ g_a(X,Y) \... 1 Provided you have raw data you could use in this process, one could use the various different models and treat them as basis functions of sorts that you wish to merge together in a least square sense. You could then merge the various models using a least square fit based on whatever data you have at your disposal. This is certainly different than simple ... 1 I guess that you are not actually interested in the variance, but in a confidence interval for your observable$\theta$. It should be noted that computing the confidence interval from the variance (i.e.$\hat{\theta}\pm 2\sigma(\theta)$) is only guaranteed to work when you estimate$\theta$with a Maximum-Likelihood estimator. The problem with your "block" ... 1 Per edited question: this is super simple. One way to show it is:$p(a,c|b) = \frac{p(a,b,c)}{p(b)} = \frac{p(a,b)}{p(b)}p(c|b) = p(a|b) p(c|b)\$