title: How to compute the standard error of the mean of an AR(1) process?

text: I try to compute the standard error of the mean for a demeaned AR(1) process $x_{t+1} = \rho x_t + \varepsilon_{t+1} =\sum\limits_{i=0}^{\infty} \rho^i \varepsilon_{t+1-i}$ Here is what I did: $$ \begin{align*} Var(\overline{x}) &= Var\left(\frac{1}{N} \sum\limits_{t=0}^{N-1} x_t\right) \\\ &= Var\left(\frac{1}{N} \sum\limits_{t=0}^{N-1} \sum\limits_{i=0}^{\infty} \rho^i \varepsilon_{t-i}\right) \\\ &= \frac{1}{N^2} Var\begin{pmatrix} \rho^0 \varepsilon_0 + & \rho^1 \varepsilon_{-1} + & \rho^2 \varepsilon_{-2} + & \cdots & \rho^{\infty} \varepsilon_{-\infty} + \\\ \rho^0 \varepsilon_1 + & \rho^1 \varepsilon_{0} + & \rho^2 \varepsilon_{-1} + & \cdots & \rho^{\infty} \varepsilon_{1-\infty} + \\\ \vdots & \vdots & \vdots & \ddots & \vdots \\\ \rho^0\varepsilon_{N-1} + & \rho^1 \varepsilon_{N-2} + & \rho^2 \varepsilon_{N-3} + & \cdots & \rho^{\infty} \varepsilon_{N-1-\infty} + \\\ \end{pmatrix} \\\ &= \frac{1}{N^2} Var\begin{pmatrix} \rho^0 \varepsilon_{N-1} + \\\ (\rho^0 + \rho^1) \varepsilon_{N-2} + \\\ (\rho^0 + \rho^1 + \rho^2) \varepsilon_{N-3} + \\\ \cdots \\\ (\rho^0 + \rho^1 + \rho^2 + \dots + \rho^{N-2}) \varepsilon_{1} + \\\ (\rho^0 + \rho^1 + \rho^2 + \dots + \rho^{N-1}) \varepsilon_{0} + \\\ (\rho^1 + \rho^2 + \rho^3 + \dots + \rho^{N}) \varepsilon_{-1} + \\\ (\rho^2 + \rho^3 + \rho^4 + \dots + \rho^{N+1}) \varepsilon_{-2} + \\\ \cdots\\\ \end{pmatrix} \\\ &= \frac{\sigma_{\varepsilon}^2}{N^2} \begin{pmatrix} \rho^0 + \\\ (\rho^0 + \rho^1) + \\\ (\rho^0 + \rho^1 + \rho^2) + \\\ \cdots \\\ (\rho^0 + \rho^1 + \rho^2 + \dots + \rho^{N-2}) + \\\ (\rho^0 + \rho^1 + \rho^2 + \dots + \rho^{N-1}) + \\\ (\rho^1 + \rho^2 + \rho^3 + \dots + \rho^{N}) + \\\ (\rho^2 + \rho^3 + \rho^4 + \dots + \rho^{N+1}) + \\\ \cdots\\\ \end{pmatrix} \\\ &= \frac{N \sigma_{\varepsilon}^2}{N^2} (\rho^0 + \rho^1 + \dots + \rho^{\infty}) \\\ &= \frac{\sigma_{\varepsilon}^2}{N} \frac{1}{1 - \rho} \\\ \end{align*} $$ Probably, not every step is done in the most obvious way, so let me add some thoughts. In the third row, I just write out to two sum-signs. Here, the matrix has N rows. In the fourth row, I realign the matrix so that there is one row for every epsilon, so the number of rows is infinite here. Note that the last three parts in the matrix have the same number of elements, just differencing by a factor $\rho$ in each row. In the fifth row, I apply the rule that the variance of the sum of independent shocks is the sum of the variances of those shocks and notice that each $\rho^j$ element is summed up $N$ times. The end result looks neat, but is probably wrong. Why do I think so? Because I run a MCS in R and things don't add up:               nrMCS <- 10000     N <- 100     pers <- 0.9     means <- numeric(nrMCS)     for (i in 1:nrMCS) {        means[i] <- mean(arima.sim(list(order=c(1,0,0), ar=pers), n = N))     }     #quantile(means, probs=c(0.025, 0.05, 0.5, 0.95, 0.975))     #That is the empirical standard error     sd(means)     0.9459876     #This should be the standard error according to my formula     1/(N*(1-pers))     0.1      Any hints on what I am doing wrong would be great! Or maybe a hint where I can find the correct derivation (I couldn't find anything). Is the problem maybe that I assume independence between the same errors? $$Var(X + X) = Var(2X) = 4Var(X) \neq 2Var(X)$$ I thought about that, but don't see where I make that erroneous assumption in my derivation. **UPDATE** I forgot to square the rhos, as Nuzhi correctly pointed out. Hence it should look like: $$ Var(\overline{x}) = \frac{\sigma_{\varepsilon}^2}{N^2} \begin{pmatrix} \rho^{2\times0} + \\\ (\rho^0 + \rho^1)^2 + \\\ (\rho^0 + \rho^1 + \rho^2)^2 + \\\ \cdots \\\ (\rho^0 + \rho^1 + \rho^2 + \dots + \rho^{N-2})^2 + \\\ (\rho^0 + \rho^1 + \rho^2 + \dots + \rho^{N-1})^2 + \\\ (\rho^1 + \rho^2 + \rho^3 + \dots + \rho^{N})^2 + \\\ (\rho^2 + \rho^3 + \rho^4 + \dots + \rho^{N+1})^2 + \\\ \cdots\\\ \end{pmatrix} $$


title: MATLAB quantiz - mean square distortion

text: have a question regarding mean square distortion returned by MATLAB's quantiz function. Here is my code:               dataset = rand(1,1000)*2-1;      M1 = 128;      qlevels1=-1+(2/M1):2/M1:1-(2/M1);      codebook1 =-1+(1/M1/2):2/M1:1-(1/M1/2);      [index1, quants1, dist1] = quantiz(dataset, qlevels1, codebook1);      Qnoise1 = sum((dataset-quants1).^2);      My after running my program my $Qnoise1$ and $dist1$ differ by magnitude of `1000` > $Qnoise1 = 0.0374$   >  $dist1 = 3.7417e-05$ Am I missing something here? Should I multiply qlevels1 by length(dataset) to get correct answer? Thanks


title: Does it make sense to average q-values?

text: I have 3 lists of significant hits (p-value < 0.05). I also adjusted for multiple comparison problem and added q-values to each list. These 3 lists (hits) also partially overlap. For the hits that overlap I was wondering if it makes sense to average the q-value? So, for instance, Hit X (List1): pval=0.004, qval=0.0045 Hit X (List2): pval=0.0006, qval=0.0006 Hit X (List3): pval=0.02, qval=0.04 Therefore, to summarize my findings I would add to "Hit X" an average q-val of (0.0045+0.0006+0.04)/3 = 0.0150. Not sure if that makes sense though?


title: Retrospective power analysis of samples from Poisson distributions

text: I understand how the simulation at Power calculation for likelihood ratio test can compute the alpha, using prop.test, and the power from a direct count of simulation values, for two Poisson distributed variables. I am interested in doing power analysis to determine the necessary number of samples, similar to using pwr.t.test(d = d, sig.level = 0.05, power = 0.8), except doing this between two samples from (suspected) Poisson distributions, and so thus not using t-tests. With given distributions I could calculate alpha and the Power here, how would I determine a suitable n? I suppose one way is to write a loop that calculates power for increasing n until the power falls within the desired threshold, but that seems computationally intense. What is the best way to approach this problem? Thanks in advance, Matt.


title: Clustering in Cox proportional hazards model MLM vs. sandwich estimator

text: This question is about a paper I am reviewing, so I cannot give a lot of detail, but I can say it involves patients clustered in hospitals and a Cox proportional hazards model. My instinct for such data would be to use multilevel modeling, since hospitals may vary on ways that are not accounted for. The authors used a robust sandwich estimator with clustering for hospitals. Is this a legitimate choice? What are its drawbacks and advantages vs. a multilevel model?


title: Can you do statistics with 4 data points?

text: I was just involved in a Q/A were a poor guy was requested to run a statistical test to _prove_ that algorithm A is better than two other algorithms. However, he has only 4 data points. Does it really make sense to make a statistical test on 4 points? Where is the limit? On three? To clarify, I understand that there are 12 numbers reported but for me it looks more like either 4 three-dimensional data points, or 3 four-dimensional data points. In their answers, authors introduce some assumptions about the underline distributions in order to artificially increase the number of data points and in the process compute mean values of four numbers; or perform t-tests on pairs of algorithms (comparing 8 numbers in total for each pair) and again making unfounded assumptions about underline distribution. How reliable is this process when you don't know the underline distribution and because you don't have enough data you can't hope to be able to infer/validate it? Isn't it more fair to just say that there isn't much you can do with so few data?


title: How do you interpret the results from ridge regression?

text: I started learning ridge regression in R. I applied the linear ridge regression to my full data set and got the following results.               gridge<-lm.ridge(divorce ~., data=divusa, lambda=seq(0,35,0.02))      select(gridge)      modified HKB estimator is 0.07693804      modified L-W estimator is 0.3088377      smallest value of GCV at 0.02 which.min(gridge$GCV)      0.02      2           round(coef(gridge)[2,-1],3)      year   unemployed femlab marriage birth military     -0.195  -0.053    0.790    0.148 -0.118   -0.042                round(coef(g)[-1],3)      year unemployed femlab marriage birth  military     -0.203  -0.049   0.808    0.150 -0.117 -0.043       **Questions:**   1. How do I interpret the results?    2. Do I have to do anything else for interpretation?


title: CDF of a random vector

text: I am reading a book that in one page it talks about cdf of a random vector. This is from the book: > Given $X=(X_1,...,X_n)$, each of the random variables $X_1, ... ,X_n$ can be > characterized from a probabilistic point of view by its cdf. However the cdf of each coordinate of a random vector does not completely describe the probabilistic behaviour of the whole vector. For instance, if $U_1$ AND $U_2$ are two independent random variables with the same cdf $G(x)$, the vectors $X=(X_1, X_2)$ defined respectively by $X_1=U_1$, $X_2=U_2$ and $X_1=U_1$, $X_2=U_1$ have each of their coordinates with the same cdf, and they are quite different. My question is: From the very last paragraph, it says $U_1$ and $U_2$ are coming from the same c.d.f. And then they define $X=(X_1, X_2)$, but they say $X=(X_1, X_2)$ is different from $X=(X_1, X_1)$. I don't really understand why the two $X$ are different. (i.e. I don't understand why $X=(X_1, X_2)$ and $X=(X_1, X_1)$ are different). Isn't $X_1$ the same as $X_2$, so it doesn't matter whether you put two $X_1$ to form $X=(X_1, X_1)$ or put one $X_1$ and one $X_2$ to form $X=(X_1, X_2)$. Shouldn't they be the same? why does the author says they are "quite different"? Could someone explain why they are different?


title: Sparse representations for denoising problems

text: I have read in a huge number of papers that sparse models (sparse coding, dictionary learning, sparse matrix factorization, ...) are good solutions for image denoising problems. I know that representing data as sparse combinations of atoms from a (overcomplete) dictionary should be the way the mammal primary visual cortex works. However it is not clear why sparse representations should eliminate noise, blur or other similar artifacts. Is there any valid mathematical explanation for this? Thinking about Lena image (384x384): ![Lena image](http://i.stack.imgur.com/EV2zA.jpg) Does its histogram really show us a "sparse nature"? ![Lena histogram](http://i.stack.imgur.com/an1vB.jpg)


title: Post-hoc comparison in glmm with interaction

text: I have data about seed predation (SP) in fruits of three differents colors (yellow, motted, dark) and in two fruiting seasons (2007, 2008). I performed a GLMM and the outcome showed that the interaction term (color:season) was significant, and some combinations of this interaction have significant Pr(>|z|), but I don't think they are the right significant combinations, because when I look the bwplot, this combinations seems to be very different from the other ones. So, I would like to know if there is any test "a posteriori" to know the p-values for each combination of color:season, and thereby be able to Know what conbination/s is/are really significant.               m1=lmer(SP ~ color + season:color +(1|Site:tree), data=datosfl, family="poisson")     AIC   BIC logLik deviance     178.3 196.6 -81.14    162.3     Random effects:     Groups      Name        Variance Std.Dev.     obsBR       (Intercept) 0.064324 0.25362      Site:tree   (Intercept) 0.266490 0.51623      Number of obs: 73, groups: obsBR, 73; Sitio:árbol, 37                              Estimate Std. Error z value Pr(>|z|)         (Intercept)            2.5089     0.2750   9.125   <2e-16 ***     colorM                -0.1140     0.3242  -0.352   0.7250         colorO                -0.6450     0.4178  -1.544   0.1227         Temporada2008         -0.7343     0.3104  -2.365   0.0180 *       colorM:Temporada2008   0.2505     0.4352   0.576   0.5648         colorO:Temporada2008   1.1445     0.5747   1.992   0.0464 *


title: Importance of variables in logistic regression

text: I am probably dealing with a problem that has probably been solved a hundred times before, but I'm not sure where to find the answer. When using logistic regression, given many features $x_1,...,x_n$ and trying to predict a binary categorical value $y$, I am interested in selecting a subset of the features which predicts $y$ well. Is there a procedure similar to the lasso that can be used? (I have only seen the lasso used for linear regression.) Is looking at the coefficients of the fitted model indicative of the importance of the different features? # Edit - Clarifications After Seeing Some of the Answers:   1. When I refer to the magnitude of the fitted coefficients, I mean those which are fitted to normalized (mean 0 and variance 1) features. Otherwise, as @probabilityislogic pointed out, 1000x would appear less important than x.   2. I am not interested in simply finding the best k-subset (as @Davide was offering), but rather weigh the importance of different features relative to each other. For example, one feature might be "age", and the other feature "age>30". Their incremental importance might be little, but both may be important. Thank you.


title: t-distribution parameter estimation

text: I know there are already several threads on this, but none seem to explicitly cover what I want. I have a set of financial data (pulled straight from Bloomberg) and am trying to fit a t-distribution (with theoretical discussion of MLE). I know it can't be solved in closed form and am looking at EM and Newton-Raphson currently, but some issues have come up:   1. I'm not sure there's any way I can justify that my data is incomplete (data is closing price of stock index over last year)   2. DF is as unknown as $\mu$ and $\sigma$, so all 3 parameters need to be estimated.   3. Data is univariate The above methods seem to need known _df_ or incomplete data. Is this strictly true? Do I need to modify the methods (if so please show) or choose new ones (if so please specify)?


title: K means clustering inadequate in determining extreme regions in R

text: ![Consider the following plot.](http://i.stack.imgur.com/QbZ5T.png) I want to identify the regions that are considerably higher than the highest cluster. (The obvious regions which should be identified as their own clusters, notably at the x coordinate ~10 e+07. How would I be able to identify that using a clustering algorithm? I am using R algorithm kmeans in the picture above. with 6 centers: `kmeans(numbers_vector, centers=6, nstart=10)` What can I do to alleviate the inadequacy of this algorithm? Use a different clustering algorithm? Have more centers? But If I have more centers, it identifies many more regions in the center (namely clusters 3,5, and 6). Any ideas? Here is histogram and density. It is important to note that the histogram and density plot DO NOT show the spike ~x=10e+07, because the number of points involved in the spike ~20, perhaps are completely overshadowed by the ~107,350 points plotted. ![enter image description here](http://i.stack.imgur.com/1zrGE.png) ![enter image description here](http://i.stack.imgur.com/Od8v7.png) Data: x (88289 obs.); Bandwidth 'bw' = 0.7574                  x                 y                  Min. : -1.272 Min. :0.000e+00   1st Qu.: 40.114 1st Qu.:4.110e-06   Median : 81.500 Median :3.167e-05   Mean : 81.500 Mean :6.035e-03   3rd Qu.:122.886 3rd Qu.:2.886e-03   Max. :164.272 Max. :4.827e-02


title: Why is the ratio of $B$ and $e^{B}$ for some variables very large?

text: I want to know how variables affect travel mode for different trip purposes (i.e. leisure trips, work trips and shopping trips) in a specified region. I have 450 respondents in three different neighborhoods. My dependent variable is mode of travel (in 5 categories including walking, public transport, taxi, private car and other modes), and also I have several independent variables such as gender, age, type of household, having car, belonging to which traffic boundaries, etc. I stratified my data by the trip purposes variable, then I used multinominal logistic regression, and my reference category was car use. However, I faced some problems in my analysis: the $B$ and $e^{B}$ in my model are very high (there are large numbers such as 43839.134 or 213.65), also the standard error is high. I don't know why, do you think it is common or not? I hope my explanation is enough.


title: To know clearly about the population

text: What is the difference among target population,study population and theoretical population ?


title: How do you report a Mann–Whitney test?

text: I am doing my dissertation, and I am conducting a number of tests. After using a Kruskal–Wallis test, I usually report the result like this: > There is a significant difference $(\chi^2_{(2)}=7.448, p=.024)$ between the > means of... But now I conducted a Mann–Whitney test, and I am not sure which values to present. SPSS gives me a Mann–Whitney $U$, Wilcoxon $W$, $Z$ and $P$-value. Do I present all these 4 values? Or are some irrelevant?


title: Minimizing the norm of a vector of parameters

text: I'm reading a paper that defines a function $f_w(x)$ that takes input $x$ and parameters $w$ and a set of constraints. There are also training data. The aim is to find the set of parameters $w$ that respect the constraints. The authors then defined the optimization problem as the following: $min_w F(w) = ||w||^2$ with respect to the constraints. My question is thus, why minimizing the squared norm of the parameters vector $w$ returns the optimal solution? Why minimization? why the norm of the parameters vector?


title: Example how maximizing and minimizing a function can be equivalent?

text: I don't understand how sometimes given an optimization problem, a function could get its optimal solution by minimizing or sometimes just by reformulation it becomes maximizing. Can you please give me an example for this?


title: Requirements to perform a two-way factorial ANCOVA

text: This is a question that landed on my desk, and I don't have the requisite experience to reply. A researcher has been asked to perform (by a journal reviewer) a two-way factorial ANCOVA on some microrarray data (a single array will be a single measurement of ~30,000 entities from a single sample in this case). The suggestion was to use ethnicity and age as co-variables. I understand that ANCOVA requires at least one categorical and one continuous predictor variables. The researcher has 3 ethnic groupings in his samples - 6 of group A, 2 of group B and 1 of group C. These samples are spread across three treatments (X, Y and Z), with 3 samples per treatment. This means the treatments cannot be split evenly by ethnicity, group C will only appear in _one_ of the treatments. Is this sample size/grouping issue a barrier to running a two-way factorial ANCOVA? Is there another test that might be more appropriate under the circumstances?


title: Simple linear model with autocorrelated errors in R

text: How do I fit a linear model with autocorrelated errors in R? In stata I would use the `prais` command, but I can't find an R equivalent...


title: How to calculate a confidence level for a Poisson distribution?

text: Would like to know how confident I can be in my $\lambda$. Anyone know of a way to set upper and lower confidence levels for a Poisson distribution?   * Observations ($n$) = 88    * Sample mean ($\lambda$) = 47.18182 what would the 95% confidence look like for this?


title: Validating a predicted standard deviation

text: I would like to validate the standard deviation which a model reports, because I think it is underestimated, because in the procedure some model selection is involved. After running the model on a validation set, it became clear to me that the standard deviation were indeed grossly underestimated. I could not really find a good way to prove this, but I think the clearest way to illustrate is a qqplot, which I generated using R:               qqnorm((m-v)/sd,ylim=c(-5,5))      with   * v (validation set)    * m (model prediction)    * sd (predicted standard deviation) And which looks like (yes, I know it is really bad): ![qqplot](http://i.stack.imgur.com/YqU4t.png) Although I find this graph very clear, it seems I was the only one in the lab, so I was wondering if anyone has suggestions on how I could plot/prove this more clearly. Also references on the topic are welcome, I made up this derivation myself because I could not find studies where also the predicted standard deviation was checked using a validation set.


title: In which case does FCM membership converge to 1/K?

text: I have tested the fuzzy C-means (FCM) algorithm using the R function `fanny` from the `cluster` package and I have wrote my own FCM algorithm to have more control on the distance function. The problem is that my final membership function converges many times to the 1/K, I mean almost all the values in the membership matrix are near to 1/K. I am wondering what are the conditions that influence the final result. I am using fuzzy parameter m=2.


title: What software is used for maps of the US (or other arbitrary areas)?

text: What software does one use to create a map that describes conditions over an arbitrary area ? The example I am thinking of would be a map of the US where states that voted one way or another would be red or blue. Obviously I can get a map of the US and use a graphics program to manually paint in the states, but I am figuring there should be some sort of software to do that in a somewhat automated fashion. Apart from a for-pay web app, I haven't been able to find anything. I am very new to the field so I'm not familiar with the available tools.


title: Modeling response times

text: I am trying to model the time lapse between when a user sees an ad and when they call the advertiser (presuming they do). I have two issues - part of the data seems exponential, but I am wondering if there is a similar distribution with an extra parameter because I cannot quite get it to fit. Also, the peak does not occur at t = 0 but rather at t = 30s or so - obviously people do not call immediately, as there is a necessary physical delay between seeing the ad and picking up the phone. But some people do call at 2 seconds, or 3, etc, and as I mentioned, the peak is at 30. Is there a simple way to model this?


title: Testing for Spatial Autocorrelation in a Negative Binomial Regression Model

text: Currently, I am working on a count data set measuring events in a third-order administrative unit. The frequency of events varies for each third-order administrative unit to account for the variability, I am using a negative binomial regression model. A crude graphical analysis of the residual of the regression indicates that no spatial autocorrelation exists. I would like to find a method to assess spatial autocorrelation in a negative binomial regression model, since a regular test of spatial autocorrelation are not feasible (ex. Global Moran's I). Do you have any suggestions?


title: Using dmnorm in JAGS

text: I have a question regarding the use of dmnorm in Jags when none of the mean and variance covariance matrix are known and fixed. I keep getting the error "node Q invalid parent value" I am not sure where the issue comes from. Any idea will be more than welcome! Thanks!               model {         for (i in 1:n) {             mu[i] <- a+b*X[i]         }         va <- (c*VM+diag)         va.p <- inverse(va)         Q ~ dmnorm(mu[], va.p[ , ])         a ~ dnorm(0,0.01)         b ~ dnorm(0,0.01)                       c ~ dunif(0.0001,100)     }


title: Simulating non-causal time series?

text: I'm looking into the possibility of using a non-causal time series filter for some data. The goal is filtering (for the purpose of anomaly detection). However, this is not particularly relevant. I'm wondering how simulation of such a time series would work? Say for example your time series model were: $$y_t = y_{t-1}+y_{t+1}+\eta_t$$ where $\eta$ corresponds to some standard random noise ($N(0,1)$ or something). My initial thoughts are to just simulate the random noise $\eta_t^*$ and then solve a linear system to get the $y_t$. However, this also sounds slow.. and I can't help but feel someone has probably already thought about this problem. So, what is the standard practice here?


title: Showing that 0 is positive recurrent

text: **The Problem** Consider the aging chain on $\\{0, 1, 2, \dots\\}$ in which for any $n \geq 0$ the individual gets 1 day older from $n$ to $n+1$ with probability $p_n$ but dies and returns to age 0 with probability $1 - p_n$. Find conditions on the $p_n$'s that guarantee that: (a) 0 is recurrent **My Work** We'll consider first the finite state space $\mathcal{S} = \\{0,1,2,\dots,N\\}$. In this finite state space, we have irreducibility and aperiodicity, so $\rho_{x0} + \rho_{xN} > 0$ for all $ x \in C = \cal S - \\{0,N\\}$. So $p(i, i+1) = p_i, p(i, 0) = 1 - p_i$. And for $h(0) = 0, h(N) = 1$, we have \begin{align*} h_N(N - 1) &= p_{N-1} \\\ h_N(N-1) &= p_{N-2}p_{N-1} \\\ & \vdots \\\ h_N(1) &= p_1p_2\cdots p_{N-1} = {\bf P}_1(T_N < T_1) \\\ {\bf P}_0(T_0 < T_N) &= p_0(1 - h_N(1)) + (1 - p_0) \\\ \xrightarrow{N \to \infty} {\bf P}_0(T_0 < \infty) &= 1 \iff \lim_{N \to \infty}h_N(1) = 0 \end{align*} Thus in order for $ {\bf P}_x(T_0 < \infty) = 1 $, we require that $\prod_{i=1}^\infty p_i = 0$, in words: that the infinite product goes to 0. (b) 0 is positive recurrent **My Work** In order for 0 to be positive recurrent, we need ${\bf E}_0T_0 < \infty$, where ${\bf E}_0T_0$ is the expected time of first return to 0. So let $\cal S = \\{0,1,2, \dots, N\\}$, let $C = S - A, A = \\{0,N\\}$, and set $g(a) = 0 \ \forall a \in A$. {\color{red}So the calculations are} \begin{align*} {\bf E}_0T_0 &= \sum_{k = 0}^\infty k {\bf P}(T_0 = k) \\\ &= \sum_{k=0}^\infty k \prod_{i=0}^{k-1} p_i(1 - p_k) \end{align*} **My Question** I am stuck with showing positive recurrence. I am also unsure how to _compute ${\bf E}_0T_0$ "directly"_ as I have been advised for a problem of this sort. And even if I were show show positive recurrence (which would imply the existence of a stationary distribution $\pi$, how would I find $\pi(0)$? Thanks for any help.


title: Does variance predict skewness?

text: I am trying to get some ideas on how to test for an implicit relationship, if any, between variance and skewness. That is, given a very large data set (e.g 90 years of monthly returns), is there a way to generally test if skewness is likely to increase with increasing variance or vice versa? I will be very grateful for any ideas or suggestions regarding this topic.


title: Variance explained for binomial proportions

text: In Rosenthal and Rubin (1979) ``A Note on Percent Variance Explained as A Measure of the Importance of Effects'', they give an example of where $r^2$ is deceptively low: > Suppose half the patients in a medical study are randomly assigned to a new > medical treatment (X = 1) while the other half are assigned the standard > medical treatment (X = 0). The dependent variable is “alive one year after > treatment” (Y = 1) vs. “dead” (Y = 0). Suppose the study obtained the > frequencies displayed in Table 1. Under the standard medical treatment, > 30%of the patients live; while under the new treatment, 70% of the patients > live. This certainly is a dramatic and important effect. Yet, for these data > r2 = .16. The conclusion that the treatment is unimportant because it > accounts for only 16%of the variance is simply wrong. Percent variance > explaified can, in some cases, then, be a very deceptive measure. The numbers are                     | Dead (Y=0) | Alive (Y=1) |     -------------------------------------     X = 0 | 35         | 15          | 50        X = 1 | 15         | 35          | 50     -------------------------------------           | 50         | 50          | 100      How do you actually calculate the $r^2$ here?


title: Discrete analog of CDF: "cumulative mass function"?

text: We call the integral of a probability density function (PDF) a cumulative distribution function (CDF). But what's the cumulative sum of a probability mass function (PMF) called? I've never heard the term "cumulative mass function" before, and the Wikipedia page for it redirects to the CDF page, so I'm confused what the proper terminology is.


title: whether to rescale indicator / binary / dummy predictors for LASSO

text: For the LASSO (and other model selecting procedures) it is crucial to rescale the predictors. The general recommendation I follow is simply to use a 0 mean, 1 standard deviation normalization for continuous variables. But what is there to do with dummies? E.g. some applied examples from the same (excellent) summer school I linked to rescales continuous variables to be between 0 and 1 (not great with outliers though), probably to be comparable to the dummies. But even that does not guarantee that the coefficients should be the same order of magnitude, and thus penalized similarly, the key reason for rescaling, no?


title: Model specification

text: Actually I want to see the impact of Investors sentiment on Aggregate stock market return, so my major concern is to identify the impact of investors sentiment. But besides investors sentiment I want to include 5 more Independent variable which might have impact also stock market return. Its a single STOCK MARKET analysis (USA), here N=1 So can anybody suggests me what model I can use to test it other than Multiple Regression Analysis?


title: How can to compare 1750 samples between 3 groups by R?

text: I have 1750 proteins that I want to compare the expression level of them between 3 groups (cell-type) using R. How can I do it?


title: How to mix probability estimators of the same phenomenon?

text: Also posted here and here. I have the following problem: I have N models that give me an estimation of the probability distribution function p(x) of a certain phenomenon x. Let's call them: $p_1(x),...,p_N(x)$. They come from different sources of information, so they can give different values, but they all refere to the same observable fact. Is there a formal and valid way to combine them into a single formula? I have read about mixture distributions, are they applicable to my case? Suppose that my models can be more or less reliable, depending on the source of information they are based on, can I also combine them giving more weight to one model rather than another? One possible solution I have thought is to make a weighted average of all the PDFs: $p(x) = w_1*p_1(x) + ... + w_N * p_N(x)$ where $\sum\limits_{i=0}^N w_i = 1$, does it make sense? Thanks very much for your suggestions! PS: to be more concrete my models give me the probability distribution that a certain person is in the position (x,y) and they rely on different sources of information like the power of a received signal or some other observable fact.


title: What does "statistically insignificantly worse" mean?

text: I am reviewing a paper where I found this sentence repeatedly: > A is statistically insignificantly worse than B I consider this to be confusing since, for me at least, it is not immediately clear if it is the same as: > A is not statistically significantly worse than B Maybe the authors are trying to convey the fact that the data for A are worse, but that this is not statistically significant. If this is the case though, it is not clear whether they mean that the sample average is worse or sample median is worse or something else. Should I recommend that they reformulate this sentence or is it common practice to say "statistically insignificantly worse"?


title: What are the limitations of Kernel methods and when to use kernel methods?

text: Kernel methods are very effective in many supervised classification tasks. So what are the limitations of kernel methods and when to use kernel methods? Especially in the large scale data era, what are the advances of kernel methods? What is the difference between kernel methods and multiple instance learning? If the data is `500x10000`, `500` is the count of samples, and `10000` is the dimension of each feature, then in this circumstance, can we use the kernel methods?


title: What test is suitable to calculate the correlation between two matrices of spatially distributed variables?

text: I have two matrices which each represent the distribution of a different variable on the same spatial domain (something like the distribution of cows on a field on one day versus the distribution of horses on the same field on a different day), and I'd like to check whether these matrices are correlated. I've been searching for an appropriate test but am still not sure what to use, so any advice would be much appreciated. Thanks a lot!


title: When is multiple comparison correction necessary?

text: I have a sort of philosophical question about when multiple comparison correction is necessary. I am measuring a continuous time varying signal (at discrete time points). Seperate events take place from time to time and I would like to establish if these events have a significant effect on the measured signal. So I can take the mean signal that follows an event, and usually I can see some effect there with a certain peak. If I choose the time of that peak and do say a t-test to determine if it is significant vs when the event doesn't occur do I need to do multiple comparison correction? Although I only ever performed one t-test (calculated 1 value), in my initial visual inspection I selected for the one with the largest potential effect from the (say) 15 different post delay time points I plotted. So do I need to do multiple comparison correction for those 15 tests I never performed? If I didn't use visual inspection, but just did the test at each event lag and choose the highest one, I surely would need to correct. I am just a little confused as to whether I do need to or not if the 'best delay' selection is made by some other criterion than the test itself (eg visual selection, highest mean etc.)


title: Why is AIC from multinom and corresponding glm different?

text: I'm reading the example of the book: "S-PLUS (and R) Manual to Accompany Agresti’s Categorical Data Analysis (2002) 2nd edition " (page 55) and when I try to reproduce the example I get this:               snoring<-c(0,2,4,5)           logit.irls<-glm(cbind(yes=c(24,35,21,30), no=c(1355,603,192,224))~snoring,      family=binomial(link=logit))           summary(logit.irls)          Null deviance: 65.9045  on 3  degrees of freedom     Residual deviance:  2.8089  on 2  degrees of freedom     **AIC: 27.061**      Ok, I get a Akaike information criterion of 27.061. But when I use the multinorm function from the "nnet" library I get this:               modLogit = multinom(cbind(c(24,35,21,30),c(1355,603,192,224))~snoring)          summary(modLogit)          Coefficients:       (Intercept)    snoring     2    3.866248 -0.3973366          Std. Errors:       (Intercept)    snoring     2   0.1662144 0.05001066          Residual Deviance: 837.7316      **AIC: 841.7316**       Why the AIC (and the Residual Deviance) are too big?


title: Classification score: SVM

text: I am using libsvm (which is meant for solving binary classification problems) for multi-class classification. How can I get classification scores / confidences for each class to effectively compare them given that libsvm can only produces scores for two classes. Desired output:               Class 1: score1          Class 2: score2          Class 3: score3          Class 4: score4


title: Testing the significance of differences between ratios with small samples

text: This is for a friend of mine. As an econometrician used to rely on large samples for inference, I find myself unsure whether the answer I have in mind is the best. Suppose we have four continuous random variables, X, Y, W and Z. We have a small (5 to 10 observations) iid sample from each. We want to test whether $\frac{E(X)}{E(Y)} = \frac{E(W)}{E(Z)}$. My first reaction is to consider taking the logs, so as to think in terms of differences rather than ratios; but I fear this would lower the power of tests that rely on the normality of the variables (X, Y, W, Z are probably approximately normal; but then their logarithm is quite far from normal). What are your thoughts on this? Any advice would be greatly appreciated!


title: How to choose the right optimization algorithm?

text: I need to find the minimum of a function. Reading the docs at http://docs.scipy.org/doc/scipy/reference/optimize.html I see that there are several algorithms that do the same thing, i.e. find the minimum. How do I know which one I should choose? some of the algorithm listed   * Minimize a function using the downhill simplex algorithm.   * Minimize a function using the BFGS algorithm.   * Minimize a function with nonlinear conjugate gradient algorithm.   * Minimize the function f using the Newton-CG method.    * Minimize a function using modified Powell's method. My function is linear. dimensionality is around 232750(this is how many different gradients I have to compute each time), it takes about 2 mins to compute the gradient and the cost once, so not cheap. I don't think I have constraints. it is deterministic and continuous.


title: Inference in undirected graphical models.

text: I know Bayesian or Directed Graphical models are good for fast inference using message passing techniques. But how does it work with undirected models? With directed models, you moralize the graph making it an undirected and triangular graph whose cliques are connected in a junction tree thus allowing the message passing. If you start out with a triangulated undirected model, you just go right to the junction tree and skip the moralization. But what if it is not triangulated? If you triangulate this, it becomes a different model. Can you do fast inference with undirected models that are not triangulated?


title: Curve fitting and max slope calculation

text: I have a dataset:   * X variable is date (from April to October)   * Y variable is vegetation biomass data In my study area, growing season starts around April when vegetation biomass is low and peaks around at the end of August when biomass is highest, and finishes around October. The purpose is to determine the exactly date when the vegetation biomass increased at maximum during the start of growing season. It should be in April. First, I did the curve fitting using Sigmoid, logistic, 4 parameter method which was the best fit for this dataset. And I got a formula of sigmoid, logistic 4 parameter method. Now, how can I calculate maximum increase date of vegetation biomass from curve fitting line to accomplish the purpose as mentioned above? Thanks a lot.


title: Why is the tick marker for zero after the bar in this qplot bar chart?

text: Feature or bug? Why is it that the tick marker for zero projects is **_after_** the bar that represents the count for zero in this plot (instead of being in the middle as I'd have expected):               > qplot(projects,data=subset(df,projects<1000),geom="bar")     stat_bin: binwidth defaulted to range/30. Use 'binwidth = x' to adjust this      ![enter image description here](http://i.stack.imgur.com/o78FC.png) Here is the data I am using:                     username gender   id   tenure  projects post     1      foo     male     123   1566      120   75     2      bar     male     456   1565       78    1     3      baz     female   678   1564       55    1


title: Adjustment of medians in the Siegel-Tukey test vs. the Wilcoxon test

text: I have two (unpaired) distributions. The distributions are far from being normal, so I'm using the Wilcoxon test to compare their medians: the Wilcoxon test does not reject the null-hypothesis (p-value = 0.5584). Then I'm comparing variances of these distributions using the Siegel-Tukey test, which turns out to reject the null-hypothesis (p-value = 0.03351) suggesting difference in variances. However, since the literature says that difference in medians might have affected the results of the Siegel-Tukey test, I've decided to rerun the Siegel-Tukey test with median adjustment (despite the previous result of the Wilcoxon test). Surprisingly, after the median adjustment the Siegel-Tukey test can no longer reject the null hypothesis (p-value = 0.1144), i.e., we can no longer claim difference in medians. Does this mean that I should have relied on the Wilcoxon tests and that adjusting the means was wrong in this case? Does this mean that the Siegel- Tukey test is not powerful enough to distinguish difference in the variances in the second case? ![Violin plots and boxplots of the two distributions.](http://i.stack.imgur.com/6zWAb.jpg)


title: Financial exposures modelling with graph theory tools

text: I was wondering how finance folks go about storing and modelling portfolio exposure relationships with the aim to later aggregate or slice & dice the exposures by different factor sets. For example, a portfolio invests into changing set of instruments (stocks, bonds, other portfolios ...) with time-varying value. The instruments in turn are sensitive to factors such as sectors (e.g. GICS classification), regions/countries, size etc. ## Tabular Approach One straightforward approach would be to express this mix of network and/or hierarchical relationships in a **tabular** format such as below, which is amenable to use with SQL design patterns (`adjacency-list`, `path- enumeration`, `nested-sets`, `closure-table`) :               > HOLDINGS        Portfolio                         Holding    Instrument       Date BALANCE.USD     1        ABC           Stock 1 Share Class A       Stock 1 2013-12-31    25360291     2        ABC           Stock 1 Share Class A       Stock 1 2014-01-31    25302011     3        ABC           Stock 1 Share Class B       Stock 1 2013-12-31    12264011     4        ABC           Stock 1 Share Class B       Stock 1 2014-01-31    12893201     5        DEF Fund 1 Share Class EUR Series 1        Fund 1 2013-12-31    21012222     6        DEF Fund 1 Share Class EUR Series 1        Fund 1 2014-01-31    21632101     7        DEF Fund 1 Share Class EUR Series 2        Fund 1 2013-12-31     8214325     8        DEF Fund 1 Share Class EUR Series 2        Fund 1 2014-01-31     8292630     9        DEF           Portfolio ABC Account Portfolio ABC 2013-12-31   155364592     10       DEF           Portfolio ABC Account Portfolio ABC 2014-01-31   156202162          > FACTORS       Instrument                                         Factor ExposureStrength     1    Stock 1                              North America: US             1.00     2    Stock 1                                    Industrials             1.00     3     Fund 1                                 Liquidity: Low             0.05     4     Fund 1                                  North America             0.70     5     Fund 1                      Europe: Eurozone: Germany             0.20     6     Fund 1 Industrials : Capital Goods: Building Products             0.25      ## Network Approach Does anyone have any experience with modelling this particular domain as **network** objects rather than tables/matrices? In `R` for example using packages such as `igraph`, or `statnet` umbrella of packages (`network`, `networkDynamic` as exposures are time-varying attributes) If so, could you kindly share short snippets of code demonstrating that it makes sense to use graph theory for these kinds of problems. My web research was not fruitful in this regard; I found mostly social network models. If not, does it make any sense at all to persist this as network objects rather than plain old tables, in terms of complexity (code overhead), computational and storage efficiency? To the support of the latter, vignette for the `network` package in R states: > For example, a network with 100,000 vertices and 100,000 edges currently > consumes approximately 74MB of RAM (R 2.6.1), versus approximately 40GB for > a full sociomatrix (a savings of approximately 99.8%). When dealing with > extremely large, sparse graphs it therefore follows that network objects are > substantially more eﬃcient than simpler representations such as adjacency > matrices.


title: Measures of predictive power of attributes in data mining

text: What are the most widely used measures of predictive power of attributes in scoring models? **Motivation** : I have a lot of attributes, more than I can study by myself and I want to select somehow the most promising ones. Is IV a good criterion for that? Are there any alternatives?


title: Is Freeman — Tukey's transformation the most powerful for percentages?

text: I've described a typical design for my experiments in this question. Well, 1-way RM ANOVA assumes a Gaussian distributed vector.   I try $y=\arcsin{\sqrt{x}}$. But for some data it works, for some doesn't... So the next step is to find the most powerful/general statistics. Freeman — Tukey's transformation with $y=\sqrt{x}+\sqrt{x+1}$ is the most powerful and seems to be appropriate in almost all my cases but again it doesn't give normal data sometimes.   What should I do this case? Box — Cox transformation?   Is it OK to transform data using different statistics within 1 thesis/article?


title: R function PRCOMP Doesn't project my 2D cloud onto the principal vector as expected

text: I have generated 100 2D correlated MVN variables in R, on which I run `prcomp`. When I plot the projected points along the first principal component (in the original coordinates) with the original data overlaid, I have a bunch of points along a line (as expected) yet they do not seem to correspond to the original points in the scatterplot (I was expecting each point along the line to be the projection of some original point) Here is my code:               require(MASS)     x <- mvrnorm(100, mu = c(0,0), Sigma = matrix(c(1,-.85,-.85,5), 2, 2),                   empirical = FALSE)       pcX <- prcomp(x, retx = TRUE, scale = FALSE)       transformed <- pcX$x[,1] %*% t( pcX$rotation[,1] )       plot(transformed, col = "red")       points(x, col = "green")        Is there some scaling going on which is preventing me from recovering the original data, or is my understanding of PCA (or R) lacking?


title: Dirac measure integral

text: I'm pretty much stuck on part (ii). Any guidance would be appreciated! The Dirac measure $\delta_{x_0}$ for $x_0$∈$\mathbf{R}$ is defined over the measurable space $\left(\mathbf{R}, \mathbf{B}\right)$ by: $\begin{equation} \delta_{x_0}(B)=\begin{cases} 1, & \text{if $x\in B$}.\\\ 0, & \text{if $x\in B^C$}. \end{cases} \end{equation}$ (i) Show that $\delta_{x_0}$ is a probability measure for all $x\in \mathbf{R}$. (ii) Show that, for every continuous function $f\colon \mathbf{R}\to [0,∞)$, $\lim_{x \to \infty}\int_\mathbf{R} f d\delta_{1/n}= \int_\mathbf{R} f d\delta_0$


title: "One-tailed" Levene Test

text: F-tests can be **two-tailed** (to test that $s_1^2 \ne s_2^2$) or **one- tailed** (to test that $s_1^2 > s_2^2$). How can I modify Levene/Brown-Forsythe to be "one-tailed", that is, to test $s_1^2 > s_2^2$ instead of $s_1^2 \ne s_2^2$? Here is a demo: ![demo](http://i.stack.imgur.com/ecd7u.png) The image shows normally distributed training data (n=1000) and a model. An F-test is used to compare the variance of one point's residuals (n=2) to the variance of _all_ of the residuals (n=2000), so the point is an outlier if its residual variance is "too large." The points are colored by p-value, where light points fit the model and dark points are outliers, and you can see that the two-tailed Brown-Forsythe rejects points that are _too close_ to the model as well as too far. _Note: A different non-parametric, one-sided variance test would be fine as well._ * * * Glen_b gave the information I needed, but I thought I would leave some implementation details (using scipy).               #basic F-test     F = var(a) / var(b)     Fp = stats.f.sf(F, df1, df2)          #Brown Forsythe     BF, BFp = stats.levene(a, b, center='median')          #two tailed t-test on transformed data     za = abs(a-median(a))     zb = abs(b-median(b))     t, tp_two_tailed = stats.ttest_ind(za, zb)          #the two tailed t test recapitulates the BF test     assert(t**2 == BF)     assert(p_BF == p_two_tailed)          #one tailed t test p value     tp = stats.t.sf(t, df)      ![scatter plots](http://i.stack.imgur.com/DPTQ3.png) Above shows scatter plots of the p values from the one-tailed $F$-test and two-tailed BF-test (left), and the one-tailed $t$ tests (right). Red points are "too close" ($s_1^2 < s_2^2$).


title: Post-stratification & quantitative variables

text: I'm in charge of contacting customers of a company in order to analyse their satisfaction. The problem is I contact them by phone and the people I contact (the sample) are not representative of the full population. Then I consider post-stratification but the problem is I need to ensure the new version of sample is representative on several quantitative variables at the time :   * age   * spent amount. I know that I need a qualitative variable for the stratification. **How can I put these 2 quantitative variables into a single qualitative stratification variable ?** **1st idea :** Split the sample based on the age quartile. I obtain 4 groups.   Split each of the 4 groups based on its spent amount quartile. Now I have 16 groups (6.25% of the population each) which can be used for the stratification **2nd idea :** Perform a clustering and find k-groups which can be used for the stratification Which one is the most used when analysts need to post-stratificate a sample based on several quantitative variables ?


title: Concavity of log likelihood for hidden markov models

text: Could you give me a good link where the concept of concavity of the log likelihood related to hidden markov model EM algorithm is clarified? Thank you in advance.


title: Confidence interval in EVT

text: I'm working through a Power-Point presentation about extreme value theory with application to finance. My question is about a technique to calculate the confidence interval of a $k$ $n$-block return level $R_{n,k}$. Suppose we have divided our data in $m$ blocks, each block of size $n$. We denote with $M^i_n$ the maximum of the $i$-th block. The Fisher-Tippett theorem tells us that we should model the distribution function of such a $M^i_n$ with a GEV distribution $H_{\xi,\mu,\sigma}(x)$. Suppose this is done using the likelihood approach, i.e. we get $H_{\hat{\xi},\hat{\mu},\hat{\sigma}}(x)$. The $R_{n,k}$ is defined to be $P(M_n>R_{n,k})=\frac{1}{k}$. As an approximation we use: $$R_{n,k}\approx H^{-1}_{\hat{\xi},\hat{\mu},\hat{\sigma}}(1-\frac{1}{k})$$ Now in one example we use some data of the S&P 500 with $n=260$ and $k=40$ with the result $$R_{260,40}\approx H^{-1}_{\hat{\xi},\hat{\mu},\hat{\sigma}}(1-\frac{1}{40})=6.8$$ which is a estimate of $R_{260,40}$. On the slides there is a sentence, I quote: "It is important to construct confidence intervals for such statistics. We use asymptotic likelihood ratio ideas to construct asymmetric intervals – i.e. the profile likelihood method." How exactly is this confidence interval for $R_{260,40}$ obtained?


title: Comparing two probabilities from the same normal distribution

text: have a normal distribution. I would like to compare two input probabilities from this population to measure how "similar" they are. Everything is subjective but I wanted to be able to say that $x$ is more "similar" to $y$ than $z$ to $y$, using some sort of an equation against the normal distribution. For example, if the mean of my population is 10, and my standard deviation is 3. I would like my simple algorithm to say find that two points (19 and 17) are more similar than two other points (9 and 10) simply because it's a lot less likely to get point 17 (since it's more than two sigmas away from the mean), thus getting another random point to be near that first point with the low probability, shows higher similarity than comparing two points that occur with much equally higher probabilities. Using something like $P(X < p_1) - P(X < p_2)$ is not good enough because I may get $0$ if both points are the same. However, obtaining two points of $9$ and $9$, should score less similarity than two points ($20$ and $20$) since $20$ is a lot less likely to occur than $9$. I feel like I need to use the difference above but somehow also use the mean and sigma to formulate a similarity "formula". Is there an existing stat test that captures what I'm trying to do above? If not, does anyone have a suggestion as how I would go about solving the problem above?


title: Can a repeated measures-design be non-temporal in nature?

text: Whenever I read about repeated measures or within-subjects design in different books and on webpages, the example that always is brought up is some kind of longitudinal or repeated measures design (e.g. people are measured on the same scale several times during the day). For my experiment, I have showed people a number of different pictures in random order, and then grouped these pictures into different categories. Is this an example of a repeated measures or within-subjects design, that is, can I use a repeated measures (within- subjects) ANOVA on my data (given that it satisfies all the basic assumptions)?


title: Statistical Significance Dependent Populations

text: I am hoping to understand best way to test statistical significance between 2 dependent population groups. For example, consider a usability test. When 100 subjects were tested, 50 of them clicked (=50% click rate). However, 50 of the subjects were male, 40 of whom clicked for an 80% click rate for males. The question is that 80% statistically significant? In other words, do men click more than the population as a whole? I think I need to use a paired $t$-test, however unsure as what I would use as the mean, since these are all population proportions.


title: Generate correlated multivariate normal samples with copula

text: I've seen examples of constructing multivariate distribution with univariate marginals coupled together via a normal copula (see  Mvdc  function from `copula` package in R). However I was wondering whether I could simulate density with multivariate marginals via a normal copula? Thank you! * * * Just to clarify how I attempted to solve this problem without copula: Suppose $vec(X_i^T) \sim N_{n \times p} (0, \Sigma_i \otimes I_p ), \, i = 1, \dots, 3$, and $\Psi ( 3x3)$ is the design covariance matrix for the samples. Let L be a $3\times 3$ matrix such that $L L^T = \Psi$, then $Y = X L $, where $X = (vec(X_1^T) \cdots vec(X_3^T))$ would be the desired correlated sample. However my concern is that this implies $Y$ follows $N_{np,3} (0, V, \Psi)$, where $V$ is the common row covariance matrix. This contradicts with my original setup. Am I wrong about this?


title: Conditional probability and Bayes theorem

text: Precision Tool Company owns a five-year-old truck. After careful consideration, management has decided that there is a one in five chance that the truck will have to have major repairs within the next year. The odds are 3 to 2 that it will require minor repair. If major repairs are necessary, there is a 0.75 probability that the company will find purchase of a new truck more economical than repairing the present truck. If minor repairs are indicated, the odds are only 1 to 4 that the purchase of a new truck will prove to be justified. If no repairs are needed, Precision will keep the present truck at least through that year.   1. Draw the probability tree and use it to determine all joint probabilities.   2. What is the probability that a new truck will not be purchased within the next year?


title: libLinear matlab impementation returns abnormal results

text: Well hello, this is my first question in CrossValidated. I hope I am in the right place. I am trying to make libLinear work. I am using the interface of matlab but it does not functions normally. To be more specific I have a dataset with multiple labels and so I tries to solved it by constructing multiple binary problems. The problem I am facing has to do with the solution is not correspond to a binary problem sometimes? I experimented with random toy samples with the following code:               S = 1;     C = 1;     e = 0.1;     B = 1;     parameters = sprintf('-s %d -c %1.5f -e %1.5f -B %1.3f', S, C, e, B);          N = 2000;     X1 = rand(N,10);     X2 = rand(N*1.2,10);     X = [X1; X2];     Xlabels = [ones(N,1); zeros(N*1.2,1)];     model = train(sparse(Xlabels), sparse(X), parameters);      A random output is               model =               Parameters: 1           nr_class: 3         nr_feature: 10               bias: 1              Label: [3x1 double]                  w: [3x11 double]      How come **w** has 3 lines? It's like its a three class problem! I checked the labels returned by the model and it adds one more like -2.1475e+09 besides the normal {0,1}! I am certainly missing something here and I cannot understand why. The problem does not occur always though. Sometimes it just returns the usual **w** with dimensions 1x11. Also I skipped the parameters argument (which is the default by the way anyway) and nothing changes. For random 60 run the number of classes varied from 2 to 8 (all values included!). Any help will be appreciated


title: What does a MANOVA tell me about my data?

text: I have a dependent variable which is measured in two different areas (DV1,DV2). I want to see if the relationship my independent variable (IV) has with each DV is transferable to the other DV i.e. if we know the relationship between DV1 and IV, can this be used to successfully predict DV2, and vice versa. People have advised me to carry out a MANOVA in order to do this, with Y containing DV1 and DV2 and x containing IV. I am a little confused with MANOVA however, as I am not sure if it is telling me what I want to know. Can someone please explain exactly what a MANOVA analysis can tell me about my data? I have posted my two outputs below. Some code I have used:               Y <- cbind(DV1,DV2)     fit <- manova(Y ~ E)          summary(fit)      gives me the output:                                  Df  Pillai approx F num Df den Df    Pr(>F)         E                  1 0.40252    13.81      2     41 2.598e-05 ***     Residuals          42                                                  ---     Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1       `summary.aov(fit)` gives me:               Response 1 :                        Df Sum Sq Mean Sq F value    Pr(>F)         E                  1  1.123 1.12299  27.826 4.341e-06 ***     Residuals          42  1.695 0.04036                           ---     Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1           Response 2 :                        Df  Sum Sq Mean Sq F value   Pr(>F)        E                  1 0.41248 0.41248  7.5361 0.008862 **     Residuals          42 2.29884 0.05473                         ---     Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1       For reference, I am actually using a pgls based MANOVA using code from "Outomuro, D., D.C. Adams, and F. Johansson. 2013. Evolution of wing shape in ornamented-winged damselflies. Evolutionary Biology." But I thought for illustratve purposes it would be easier to stick to the normal non- phylogenetic MANOVA.


title: Non-linear SVM classification with RBF kernel

text: I'm implementing a non-linear SVM classifier with RBF kernel. I was told that the only difference from a normal SVM was that I had to simply replace the dot product with a kernel function: $$ K(x_i,x_j)=\exp\left(-\frac{||x_i- x_j||^2}{2\sigma^2}\right) $$ I know how a normal linear SVM works, that is, after solving the quadratic optimization problem (dual task), I compute the optimal dividing hyperplane as $$ w^*=\sum_{i \in SV} h_i y_i x_i $$ and the offset of the hyperplane $$ b^*=\frac{1}{|SV|}\sum_{i \in SV}\left(y_i - \sum_{j=1}^N\left(h_j y_j x_j^T x_i\right)\right) $$ respectively, where $x$ is a list of my training vectors, $y$ are their respective labels ($y_i \in \\{-1,1\\}$), $h$ are the Lagrangian coefficients and $SV$ is a set of support vectors. After that, I can use $w^*$ and $b^*$ alone to easily classify: $c_x=\text{sign}(w^Tx+b)$. However, I don't think I can do such a thing with an RBF kernel. I found some materials suggesting that $K(x,y)=\phi(x)\phi(y)$. That would make it easy. Nevertheless, I don't think such a decomposition exists for this kernel and it's not mentioned anywhere. Is the situation so that all the support vectors are needed for the classification? If so, how do I classify in that case?


title: Is a Stationary VAR Process with Zero Mean Gaussian Innovations a Gaussian Stationary Process?

text: Consider the stationary VAR process $${\bf X}_t = \sum_{\tau = 1}^{L} A_\tau {\bf X}_{t-\tau} +{\bf \epsilon}_t$$ If the innovations $\epsilon_t \sim MVN({\bf 0},\Sigma)$ then is ${\bf X}_t$ a Gaussian stationary process? Is it correct that due to the invertibility of the VAR into an MA and observing that ${\bf X}_t$ is the sum of zero mean MVN random variables the above is true, or is there a flaw somewhere in this.


title: In which Data Stream Mining Algorithms do Damped Windows make sense?

text: For Data Stream Mining, especially in Document Classification, the most common ML algorithms are Multinomial Naive Bayes, Stochastic Gradient Descent and Ozbag (ADWIN). When looking at their ageing behavior, they either assign same weights to all data points (Landmark Window) or "forget" older data points (Sliding Window). The idea of a Damped Window has been introduced in research, which would assume a weighting function that would make the classifier decrease the weight of data points exponentially. So far, I have not seen these damped windows implemented in any ML libraries. Is there any logical reason against using them? In which ML algorithms would they make sense?


title: Writing an equation between graph elements

text: Can we write an equation for **P(R|S) = ?** ![enter image description here](http://i.stack.imgur.com/41Dg3.png) If you refer to a document it will be great.


title: I want to show $e^{-\alpha t}B(e^{2\alpha t})$ is a Gaussian process and I find mean and covariance functions

text: Let $B(t)$ be Brownian motion. Show that $e^{-\alpha t}B(e^{2\alpha t})$ is a Gaussian process. Find its mean and covariance functions. thanks .


title: Create distribution model from samples?

text: If I have a set of samples, say 100-200 samples and I'd like to create the distribution model from this list of samples, what is a reasonably efficient way of doing it? Are there any opensource / easily accessible statistical libraries? eg: if I _assume_ it's a normal distribution, I can easily find the mean and variance and call it a day. However, the distribution might actually have non- negligible 3rd (skewness) or even 4th moment (Kurtosis) so my "assumption" is not very accurate. My gut feeling is this may work:               // assume samples already in samples[]     float avg = CalculateAverage(samples);     float variance = CalculateVariance(samples);     float skewness = CalculateSkewness(samples);     float kurtosis = CalculateKurtosis(samples);     string definedBy = "average, variance";          if(skewness > skewThreshold)        definedBy += ", skewness";     if(kurtosis > kurtosisThreshold)        definedBy += ", kurtosis";     Console.Writeline("This distribution is defined by " + definedBy);      So: Q1: Will it achieve the purpose of classifying the distribution/creating the distribution model ? Q2: What values of skewThreshold and kurtosisThreshold are reasonable?


title: Help with my time series ARX model prediction?

text: I have created an ARX-model where I predict the nitrogen oxide levels based on past values of nitrogen oxide with past exogenous input values nitrogen dioxide, temperature, atmospheric particulate matter, etc. This might be a straightforward question, but I started to wonder if for example my ARX-model would be the following: $$ Y_{t} = \alpha_{0} + \alpha_{1}Y_{t-1} + \alpha_{2}Y_{t-2} + \beta_{1}X_{t-1}^{1} + \beta_{2}X_{t-2}^{1} + \gamma_{1}X_{t-1}^{2} + \gamma_{2}X_{t-2}^{2} + . . . $$ Where Y:s correspond to my past values of nitrogen oxide(NO) and the X:s correspond to different exogenous input values. Now if I have data available until time t on all time series and I want to make predictions until say time t+10, won't I need to first predict all the future values of exogenous inputs? After I have predicted the exogenous values I can start to predict the NO-levels with my model correct? So I'm predicting predictions with predictions? :) I'm very new at this and my question might be stupid x) Please comment if my question is unclear


title: What are the use cases related to cluster analysis of different distance metrics?

text: I'm trying to use different distance metrics like Euclidean, Manhattan, cosine, chebyshev among other distance metrics in my k-means algorithm to calculate distances between the data points and the centers. In what situation would one distance metric be more useful over the other in a clustering scenario? [Comparing all the above mentioned distance metrics]


title: validating a regression model

text: I want to build a regression model to predict daily income from customers. I have 2 problems: 1\. Choosing data for the training set - do I use daily income from 1 month ago, 6 month ago etc. 2\. How do I validate the model - I can use the model to predict the income from the customer on the first day, the second day etc. On some days the model has a small prediction error and on other days large error therefor I do not know how to validate the model. Help please!


title: Packages or libraries for multiple-output learning

text: I have a learning problem from $X$ to $Y$ where:   * $X$ = $n$ input numeric vectors of $m$ dimensions    * $Y$ = $n$ output numeric vectors of $k$ dimensions In other words:       ![enter image description here](http://i.stack.imgur.com/rgcru.png) I am hoping to collect a **list** of **R** packages or **Python** libraries for **multiple-output** problems for classification and regression. For example, do any of the learning methods in caret support this functionality? What packages in general are available for this problem?


title: Multivariate linear regression model

text: We have the standard linear multivariate regression model $y=\beta_0+\beta_1x_1+\beta_2x_2+u$ under the Gauss-Markov assumptions. Suppose we estimate $\gamma_0, \gamma_1, \gamma_2$ from $x_2=\gamma_0+\gamma_1x_1+\gamma_2y+v$ and get $\hat{\gamma_0}, \hat{\gamma_1}, \hat{\gamma_2}$. Is $\frac{1}{\hat{\gamma_2}}$ an unbiased estimate of $\beta_2$? I don't know how to answer the question. First of all it seems that to run the regression of $x_2$ on $y$ and $x_1$ we need $\beta_2$ be nonzero. And is the expression $x_2=\frac{-\beta_0}{\beta_2}+\frac{-\beta_1}{\beta_2}x_1+\frac{1}{\beta_2}y+\frac{-1}{\beta_2}u$ somehow related to the regression model $x_2=\gamma_0+\gamma_1x_1+\gamma_2y+v$? Another question is whether the Gauss-Markov assumptions hold in the case of this new regression, seems like some of them can fail, for instance $E(v|x_1,y)=0$ may not be true(is it important here?). Maybe stupid questions, but I am new to this topic.


title: Repeated measure problem (Discrete variables)

text: > **Possible Duplicate:**   >  Difference between generalized linear models & generalized linear mixed > models in SPSS   >  Repeated measure problem I have 6 subjects. For each subject, we measured two variables over six-week period once a week. Let’s call the variables A and B. Variable A is binary (male or female), and variable B is also discrete (blood pressure: High, medium, and low). Now I want to see if females have more highe blood pressure levels than males. What kind of statistical method do you suggest here (keep in mind that we’re dealing with repeated measure)?


title: Difference between generalized linear models & generalized linear mixed models in SPSS

text: I am wondering what the differences are in SPSS between   * `analyze-> generalized linear models-> generalized linear models` &   * `analyze-> mixed models-> generalized linear` Do they deal with missing values differently? My dependent variable is binary and I have several categorical and continuous independent variables.


title: Time series regression - ML estimation

text: I have a linear regression model with some correlated errors: $Y_t=\beta_0+\beta_1X_1+\beta_2X_2+\epsilon_t$, where $\epsilon_t$ is a AR(1) i.e. $\epsilon_t=\phi\epsilon_{t-1}+\nu_t$ with $\nu_t$ as white noise terms. I want to obtain the joint maximum likelihood estimates of $\beta_j$'s for $j=0,1,2$ and $\phi$ as well as the residual standard error i.e. $\sigma$ for the white noise term. I guess it can be done with a time series package in R, maybe the forecast package. But I am not sure exactly how to write the code. My first guess is something like this:               library(forecast)        Arima(dat[,c("Y")],order=c(1,0,0),xreg=dat[,c("X1","X2")],           include.drift=TRUE,method=c("ML"))      Please let me know if there is any package available for this. Thanks.


title: Spearman's rank correlation for a beginner: Sample size and Interpretation

text: I am a computer scientist performing research, which includes calculating the Spearman rank correlation of two lists, one ranked by a human, another by a computer program. I have the following questions:   1. I was reading this post and have become concerned with how large my sample size has to be to have a significant confidence level.    2. Also if I have understand correctly the Fisher transform gives you a confidence level after you have gotten the Spearman rank correlation. If possible I would like to know the ideal sample size before I start the experimen   3. On interpreting the value given to you by the Spearman rank correlation: I have read about the null hypothesis > The general form of a null hypothesis for a Spearman correlation is: H0: > There is no association between the two variables [in the population]. > However what I want to prove is the opposite ie: that there is a very close association between the two variables. How is this done?


title: What are the statistic for the percentage of illegitimate baby in US?

text: Is there a way to know that? Is this a good stackexchange to ask? Basically I am looking at the number of children within marriage whose biological dad is not the mother's husband during conception.


title: Predicting Required Number of Vehicles

text: I'm trying to predict the number of vehicles per model required to maximise utilisation and booking volume but I've hit a brick wall. I've done some exploratory analysis and utilisation is affected by booking days. Also booking days is affected by booking length and bookings per car. I've worked out that 5x 4-day bookings would reach the target but the number of vehicles required depends on how popular the model is but I'm not sure where to go next.


title: Interpretation of Random Slope in Multilevel models

text: i got a question regarding random slopes in hierarchical models (multilevel). I fit two models: Model 1: without random slope for covariate x (only random intercept) Model 2: random slope for covariate x When adding the random slope (Model 2) for x the fixed effect for covariate x gets insignificant, but the fit of the model improves. How can I interpret the fixed effect? Does it imply, that the mean impact of covariate x is insignificant, when allowing different slopes for the groups? Is it possible that the impact of x is significant for some groups and insignificant for others? How can i check this? Thanks


title: How to assess whether IQ predicts emotional experiment more at time 1 or time 2 controlling for covariates?

text: I have the following data set of a two-wave longitudinal study in educational psychology. A group of 400 people did some questionnaire on IQ, personality and two different emotional states (sadness and melancholy) at time 1 and one year later at time 2. My hypotheses are:   * IQ is a significant predictor emotionalstate1 and emotionalstate2 over and above the personality factors.    * IQ is a stronger predictor for emotionalstate2 than for emotionalstate1.  My variables are:   * time 1: IQ, five personality factors, emotionalstate1, emotionalstate2    * time 2: emotionalstate1 emotionalstate2 I did multiple regression for each emotional state and each time point separately. But how could I account for the repeated measurements in my analysis? Keep in mind that the assumption of independence in GLM will not be met then.


title: Precision of replicated measurements on the same subject

text: I have developed a technique to count bacteria on a petri dish. Previously these bacteria could not be counted, so there is no reference standerd to compare our measurements with. To show the preciseness of my method I thought it was a good idea to do replicated measurements. Therefor I measured the same 200 petri dishes three times using my technique. The amount of bacteria did not change between measurements. The measurements are continues and not normally distributed. Our dataset looks like this:                    M1       M2       M3      34,543   34,500   35,001       0,002    0,002    0,002     231,003  230,974  231,509       9,456    9,693    9,543      35,543   35,500   35,001       0,034    0,035    0,035      31,603   30,674   31,609       40,456   39,693   39,543        etc       I want to know if my new method is consistent. How can I test that my three measurements of the same 200 samples do not differ significantly? If it implies a statistical test (e.g. repeated measures anova), would a p>0.05 be enough to enforce my claim?


title: manipulating Gaussian distributions

text: I have a joint distribution which factorises as follows: $$ q(w) = \prod_{i=1}^{N} t_i $$ where each $t_i$ is a 3-dimensional normal distribution which some mean $\mu_i$ and variance $\mathbf{I}\sigma_i$ where $I$ is the identity and $\sigma_i$ is the precision which is also three dimensional to have anisotropic behaviour. Now, if I want to combine this into one big multivariate distribution ($q(w)$) which is of dimension $N \times 3$ is it simply a matter of stacking up the mean and sigma vectors. I think, yes, as the factorization implies independence but wanted to check.


title: Generating a categorical variable from an imputed variable

text: I am using multiple imputation to impute a continuous variable ($X$) with $\approx30\%$ missing values. I have a question regarding the generation of a new categorical variable ($Y$), starting from this imputed variable. I want to use $Y$ as an independent variable in a logistic regression model instead of $X$ (I understand I'm losing information by doing this). This is an example of how the categories in $Y$ would work:   $$ Y = \begin{cases} 0 &{\rm if}\ \ \quad\quad\quad\\! X < 100 \\\ 1 &{\rm if}\ \ 100\le X < 200 \\\ 2 &{\rm if}\ \ 200\le X < 300 \\\ 3 &{\rm if}\ \ 300 < X \end{cases} $$ As imputed values vary between each set, I end up with a (not that much) different number of cases in each category of $Y$ on each imputed set. Is this an appropriate approach if $X$ was imputed?


title: Logistic regression with binary inputs

text: Is logistic regression an appropriate classifier when the input data are binary? Say we are conducting an experiment where the subject is presented with blue and green circles of varying shades, and asked to pick the blue one in each trial. We have 10 controls and 10 cases who are presumed to be blue- green color blind. Each subject is presented with 30 trials.We record the data in a matrix like this:               [[1 1 0 1 0 0 1 1 0 1 1 0 1 1 0 0 1 1 1 1 1 1 0 1 1 1 1 1 1 1]      [1 1 1 1 0 1 1 0 0 0 1 0 1 0 1 1 1 1 1 1 1 0 1 1 1 1 1 1 1 1]      [1 1 0 0 1 1 1 1 1 1 0 1 1 0 1 0 0 1 1 1 1 1 1 0 1 1 1 1 0 1]      [1 1 0 1 1 0 1 1 0 0 1 1 1 1 1 0 1 1 1 1 1 1 1 1 1 1 1 1 1 1]      [1 0 1 1 1 0 1 0 1 1 1 1 1 1 0 1 0 1 1 1 1 1 1 1 1 1 1 0 1 1]      [1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 0 1 0 1 1 1 1 1 1 1 1 1 1 1 1]      [1 1 1 0 1 1 1 0 1 1 1 1 1 1 1 1 1 1 1 1 1 1 0 1 1 1 1 1 1 1]      [1 1 1 1 1 0 1 0 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1]      [1 1 1 0 1 1 1 1 1 1 1 1 1 1 0 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1]      [1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1]      [0 0 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 1 0 0 0 0 0]      [0 0 0 0 0 1 0 1 1 0 0 0 0 0 0 1 0 0 0 1 0 0 0 0 1 0 0 1 0 1]      [0 0 0 0 0 0 0 1 1 0 0 0 1 0 1 0 0 0 0 1 0 1 0 0 0 0 0 1 0 1]      [1 0 1 0 1 0 1 0 0 0 0 0 0 0 0 0 0 1 0 1 0 1 0 0 1 0 0 0 0 1]      [0 1 0 1 0 0 0 0 1 1 0 0 1 1 0 0 0 0 0 1 0 0 1 1 1 0 0 0 1 1]      [0 0 1 0 0 0 1 1 0 1 1 0 1 1 0 0 0 0 0 1 0 0 1 1 0 1 0 0 1 0]      [0 0 1 0 1 1 1 0 0 0 1 0 1 1 1 1 1 0 0 0 1 1 1 0 1 1 0 0 0 1]      [0 1 0 0 1 1 1 1 0 1 0 1 0 1 1 0 0 0 0 1 1 1 1 0 1 0 0 1 1 1]      [0 0 1 0 1 0 1 0 0 1 1 1 1 1 1 1 0 1 1 0 0 1 1 1 1 1 0 0 1 1]      [1 1 1 1 1 1 0 0 0 1 1 1 1 1 1 1 1 1 1 1 1 0 1 1 1 1 1 1 1 1]]      The first 10 rows record the data for the controls, the last 10 are for the cases. From these data we can calculate mean and standard deviation of the success rate for each subject, and we can calculate the effect size between the two groups. Given that we have a "large" effect size, we would expect to see high specificity and sensitivity in a binary classifier for "color blind" vs "not color blind." Is logistic regression an appropriate classifier in this case? What should the exogenous variables be? Should I use the mean of the 1's and 0's for each subject? Should I use all 30 "features" as my exogenous variables? If the input data are in random order, i.e., exogenous variables are not actually tied to a "feature" _per se_ they are just a sequence of trials, can I still treat the trials as features and expect the logit to minimize the trial(s) that are less informative? In this answer it seems that a Chi-square test may be an alternative to logit, but it's not clear how that would work. I have implemented this classifier in Python using both statsmodel and sklearn. The statsmodel library is unable to converge when I use all 30 features, presumably because issues arise from multicollinearity. **Edit:** I wrote a naive Bayesian classifier based on this procedure, and it works quite well. It gives 75% classification accuracy, which is about what I expected based on the 1.4 effect size between my two groups. Is this the best approach, or at least a reasonable one?


title: Bivariate normal expectation of the sinus cardinal

text: I would like to get an analytical expression for $$\mathbb{E}\left(\frac{\sin(aX)}{aX}\frac{\sin(bY)}{bY}\right)$$ or at least an analytical approximation thereof, when $a,b$ are positive reals, and $$\begin{pmatrix}X\\\Y\end{pmatrix} \sim \mathcal{N}(\mu,\Sigma) $$ where $$\mu \equiv \begin{pmatrix}x^\circ \\\ y^\circ\end{pmatrix} \quad \Sigma \equiv \begin{pmatrix} 1 & \rho \\\ \rho & 1 \end{pmatrix}$$ Using either the probability density function or the characteristic function to compute the expectation leads me to an expression with a single integral, but I can't go further than this.


title: How to plot a 5D data set in "star coordinates"?

text: I am reading the paper "Star Coordinates: A Multidimensional Visualization Technique with Uniform Treatment of Dimensions" and trying to plot my data. Lets say I have $A(2,5,3,1,8)$, a five dimensional data point, and points are calculated by the formula explained in the paper. > The basic idea of Star Coordinates is to arrange the coordinate axes on a > circle on a two-dimensional plane with equal (initially) angles between the > axes with an origin at the center of the circle (Figure 1). Initially, all > axes have the same length. Data points are scaled to the length of the axis, > with the minimum mapping to the origin and the maximum to the other end of > the axis. Unit vectors are calculated accordingly. ... > > This is simply an extension of typical 2d and 3d scatter-plots to higher > dimensions with normalization. I have hard time grasping the idea. How do it plot it? The main problem is I could not understand the formula in the paper.


title: Can I use t-test for interaction in a two-way factorial ANOVA?

text: I have a two-way factorial design, with factors treatment (A vs. B) and time (T1 vs. T2). The same subject was measured twice, at T1 and T2, so I have repeated measurements. I am interested to see whether there is interaction between treatment and time on the response variable $x$. I am applying a t-test. $H_0$ is: $(\mu_{A,T2}-\mu_{A,T1})-(\mu_{B,T2}-\mu_{B,T1})=0$, that is, there is no differences in the average change (from T1 to T2) between the treatment A and B groups. To conduct the t-test, I computed the differences of $x$ between T2 and T1, this gives me $\Delta x_{A}$ and $\Delta x_{B}$, and then I can apply a t-test to see whether the average of $\Delta x$ are significantly different between the A and B groups. My question is:   1. Is this a statistically sound way to test the hypothesis?   2. What's the difference between this method and the classical ANOVA method, which computes the mean square and uses the F-test? Which one is more powerful?


title: When does drawing additional samples not cause much benefit in estimating the population probability distribution?

text: Suppose I have a random variable X that either evaluates to A, B, or C when I realize it (each with its own fixed probability). Now let's say I have drawn 100 samples and the frequency distribution is something like               A  30         B  55        C  15        dividing by n=100 gives me the probability estimates. However, when do I have a 95% chance of being $\delta$ close to the true population probabilities? EDIT: I'm not precise enough with that question. Does $\delta$ refer to the best deviation, the worst deviation or the average deviation? If some of them are tricky please just assume the simple interpretation. I just need a rule of thumb.


title: Runs Test and Chi Square Distribution

text: I want to identify random data by applying some tests to the observed byte stream. I used the chi square test already on a frequency analysis, which works fine. To reduce the false-positive rate I want to apply some more tests, like the runs test, which checks for monotonous sequences in the byte sequence (not the runs on the binary level, but on the byte level [values from 0 - 255]) I found some papers which describe the procedure, but my implementation doesn't work, probably I understood something wrong. The papers are unfortunately in German, but what they say in general, is 1) we count either the monotonically increasing or monotonically decreasing sequence of bytes (i.e. 3 | 1, 6, 9, 44, 74 | 11 | 6, 251, 46 | ... => two monotonically increasing sequences). 2) the probability, each sequence must theoretically occur, is defined by p = r / (r+1)!, where r is the length of the run. In total I have 6 categories, where category no. 6 includes all runs > 5. When I perform this test on my data, I get a chi-square distribution which is much too high, but the data is definitely random (chi square value between 80 - 100). I performed this test several times on different samples, the distribution remains constant, and the chi-square value remains constant, too, but too high. I am not 100% sure about the calculation of the probability. thank you very much for your help in advance!


title: What can be a cause of a extremely high standard coefficient?

text: I am using RapidMiner to perform linear regression with ridge parameter 1$\text{E}$ -8, min tolerance 0.05 and M5 prime for feature elimination. The `std coefficient` parameters appear to be excessively high. What are some possible causes of this?               property        user_created   user_question_count    intercept     coefficient     -0.076         0.046                  0.767     std. error      0.012          0.012                  0.004     std coefficient 5.59000E12     -9.271701831E12        ?     tolerance       0.984          0.984                  ?     t-Stat          -6.523         3.962                  199.67     p-value         0              0.000                  0


title: Robust Wald test for Poisson with Stata

text: I have a Poisson model which I use on a medical data set with 329 observations. In the regression I am particularly interested in the simple Wald statistic for a single coefficient: $$t = \frac{\widehat{\beta} - \beta}{s.e.(\widehat{\beta})}$$ I was told that simply using the robust option will not help and it was suggested that I use the percentile t-method for conducting the Wald test to refine the t-statistic. I know that I have to use the bootstrap command in Stata but I have not figured it out completely. The theoretical points about this type of refinement are clear but the problem is to implement it and then to know what the resulting $t$ tells me compared to the un-refined one. Any guidance on this would be much appreciated.


title: Different ways to calculate conditional probability of independent events

text: Say we have TWO tests for cancer. Each test has the same probabilities of being right/wrong for cancer/no cancer. What is the probability of getting two positive test results? In what follows, I'll use C for "has cancer", "+" for a positive test result, "++" for two positives (I don't distinguish between the tests because they have the same probabilities). ### A. Use the formula for total probability: P(++) = P(++|C) $\times$ P(C) $+$ P(++|$\neg$C) $\times$ P($\neg$C) ### Since each test is independent of the other, apply the formula for joint probability: = P(+|C)$^2$ $\times$ P(C) $+$ P(+|$\neg$C)$^2$ $\times$ P($\neg$C) ### B. Since each test is independent of the other one, the joint probability is the product of the probability of each (which is the same): P(++) = P(+)$^2$ ### Then, apply the formula for total probability to each: = [P(+|C) $\times$ P(C) $+$ P(+|$\neg$C) $\times$ P($\neg$C)]$^2$ Obviously, the two are not equal. Why not? And which one is correct?


title: Physical meaning of correlation?

text: Assume that X depicts a random variable denoting the time it takes someone to sweep the floor today and Y be the time it takes him tomorrow and Z be the time it takes him on the last day of October.   1. If X, Y and Z are assumed to be independent, it means that the person does not "learn" how to sweep faster through his first day's experience and the result of X is irrelevant to the result of Y and Z.    2. What does it mean physically when X,Y and Z have a correlation coeffecient of a certain amount (say 0.5). I understand that $\rho_{XY} = \dfrac{cov(X,Y)}{\sigma_X \sigma_Y} = \dfrac{E[(X-\mu_X)(Y-\mu_Y)]}{\sigma_X \sigma_Y} $


title: What reference can I cite for the proof that the sum of n exponential variables follows a gamma distribution?

text: There is a fairly common theorem, which states that: The sum of $n$ independent variables following an exponential distribution $\mathrm{Exp}(\alpha)$ follow an gamma distribution $\mathrm{Gamma} (n, 1/\alpha)$ (also known as Erlang distribution). I'm using this theorem in my thesis. I've been asked to cite the paper that introduces, or, if it's not possible, to at least cite a paper that explicitly mentions it. Anyone knows of any such papers?


title: Classifier jars for java

text: I have a train file which has categorical features like               IN  JJ  PRP_VBP     VB  NN  PRP$     . . .      The third column is the ground truth and can have value only out of {PRP_VBP, PRP$}. The first two columns can have value from 36 such tags. Now I have to train some model and then predict the ground truth in java. Which library(jar) can be used and which classifier will particularly be useful for this?


title: Optimal sampling strategy for EFA, CFA and SEM

text: I'm wondering what should be the optimal sampling strategy for my dissertation research. I have four data sources (two open source software projects meta- repositories and two global startup databases). I'd like to perform EFA to discover (or, rather, confirm my theory-based assumptions) the factor structure of the study's constructs. Then I plan to perform CFA to determine validity and reliability of the measurement model. Finally, I plan to perform SEM analysis to test the study's structural model and hypotheses. Having said that, I plan to perform data analysis (at least, SEM - not sure about EFA/CFA) on two data sets: pilot and main. I believe that pilot analysis will allow me to modify model, if fit indices will be inadequate. Then I plan to perform main SEM analysis of the modified model (and, possibly, alternative models) on the main data set. In addition, I plan to perform both covariance-based (CB) and partial least squares (PLS) SEM analysis to compare them for my study. What should be the **optimal approach and its steps** in terms of the following:   1. Sampling technique. I was thinking about randomized sampling of data from each OSS meta-repo or from a merged data set; then selecting corresponding data from startup databases.   2. Strategy on dividing the sample data set into pilot and main data sets.   3. Any special steps for sampling due to multiple methods (EFA/CFA/SEM).   4. Any special steps for sampling due to alternative models.   5. Any special steps for sampling due to analyzing via both CB-SEM and PLS-SEM. Bonus question: :-) I plan my study as cross-sectional, but the data in meta- repositories do not exist for exactly the same time frames. For example, data from seven OSS repositories are within range from September 2012 to December 2013. I think that for OSS world the variance within projects' characteristics should not be dramatic, as OSS ecosystem is not very dynamic on average. **The question is whether using this semi-cross-sectional approach will allow me to retain statistical validity and what statistical tests exist to confirm that?** You help and advice on this is greatly appreciated!


title: Estimating linear regression with OLS vs. ML

text: Assume that I'm going to estimate a linear regression where I assume $u\sim N(0,\sigma^2)$. What is the benefit of OLS against ML estimation? I know that we need to know a distribution of $u$ when we use ML Methods, but since I assume $u\sim N(0,\sigma^2)$ whether I use ML or OLS this point seems to be irrelevant. Thus the only advantage of OLS should be in the asymptotic features of the $\beta$ estimators. Or do we have other advantages of the OLS method?


title: How to calculate the similarity of two corpora (each of which contains a set of documents)?

text: I have two corpora for example, each of which contains a set of different documents, and each document are already represented as a vector of words in a certain way. The two corpora are small, only contains 20 documents or so. So How can I calculate the similarity between the two corpora? Any good ideas? I think the more overlapping words / documents between two corpora, the more similar they are. Currently I have two simple methods implemented, one is `a * DocumentVectorSimilarityOfTwocorpora + (1-a) * WordVectorSimilarityOfTwocorpora`, the other one is to represent each corpus as a document-word matrix, then unfold it into a one way vector, then use cosine similarity. Both of them are simple ways, I am wondering if there are more accurate methods. Is there a way to obtain the distributions of implicit topics of the two corpora and compute cosine similarity of the two topic distribution vectors of them?


title: Sampling from an arbitrary distribution with unknown CDF

text: I have a continuous distribution whose PDF I know the expression for but whose CDF is difficult to compute analytically. I understand that if I know the CDF value, then I can use inverse transform sampling to sample from the distribution. This would involve drawing from a uniform(0, 1) distribution, finding the corresponding point on the CDF, and inverting the CDF to get the value of the random variable with the corresponding CDF value. Moreover, the PDF is a complicated function whose shape is not always easy to visualize. If I knew the shape, I could easily set up a proposal distribution and employ rejection sampling. In this scenario, is there any way to generate samples from the distribution given only the expression for the PDF? For example, is evaluating the PDF at thousands of points along its domain, normalizing these probability values, and drawing from a multinomial distribution to get the index of the draw a valid way to draw from this distribution?


title: Feature selection for low probability event prediction

text: I'm currently trying to predict the probability for low probability events (~1%). I have large DB with ~200,000 vectors (~2000 plus examples) with ~200 features. I'm trying to find the the best features for my problem. What are the recommended method? (preferred in Python or R, but not necessarily) Thanks!


title: Marginal effects of a logistic model, and their standard errors

text: I am using SAS to estimate some logistic models. Usually, I work with either MDs or social scientists, and odds ratios are the preferred metric. But I am now working with a client in economics/law and she wants the marginal effects and their standard errors, and she wants them at the means of the other variables. This isn't easy in SAS, but, with help from tech support, I found that you can do this with PROC NLMIXED, I believe you then need the out = der option. Something like this               proc nlmixed data=olivia.small;       p=1/(1+exp(-(Intercept+ba*log_fund_age + bb*log_fund_size + bc*yield + bd*loaded + be*log_assets)));       model vote_code_num ~ binomial(1,p);       parms intercept 36.43 ba -14.55 bb -0.98 bc -0.37 bd 2.2 be -0.07;       predict p*(1-p)*ba out=a der;       predict p*(1-p)*bb out=b der;       predict p*(1-p)*bc out=c der;       predict p*(1-p)*bd out=d der;       predict p*(1-p)*be out=e der;       where year = 2003;       run;      but then the output data sets a, b, c, d and e have the derivatives and their standard deviations for each observation in the data set, not for the mean of the other variables. It's easy to find the mean of all those derivatives, but 1) Is that the same as the marginal at the mean of the other variables? and 2) How then to get the standard errors? Peter


title: Standard deviation in Weibull distribution

text: I've have a script that creates Weibull distributions with a static shape parameter and a variable scale parameter, the scale parameter changing to yield a distribution that meets a specified mean while maintaining the shape parameter. I have calculated that the standard deviation for my distributions always equals 2*sqrt(5). This is really small for what I am doing: I am using the distribution to model webpage read-time based on a known average (calculated from the number of words on the page). Times, including standard deviation, are in seconds. 2*sqrt(5) doesn't seem like a whole lot of wiggle room for read-time variation when trying to emulate an "average" web browsing user. How can I approximate how many standard deviations "wiggle room" is appropriate for my application? I was surprised, since I thought that "standard" deviation implied that a majority of samples would be within one standard deviation of the mean. I can't believe it's so tight around my average. It doesn't seem like a "standard" deviation is "standard" at all, since in other distributions I've looked at, a large percentage of samples fall within one of them. This assumes I didn't make a calculation error -- I'm not finding an error yet at least. Thanks


title: Are there any good papers comparing different philosophical views of cluster analysis?

text: Lots of people use cluster analysis. I've heard very few explicitly say why. I imagine this is because within a given field, most researchers seem to understand why clustering is used for the problems typical to that area - but uses vary between fields, and I haven't seen those differences addressed anywhere. I am particularly interested in the contrast between latent-variable interpretations of model-based methods (mixture models), and clustering applications in machine learning that don't seem as worried about interpretation of clusters, just that they are useful in some way. Data reduction is an even more agnostic application that is very common. There are loads of papers comparing different methods for clustering - but I can't find any that compare philisophical/theoretical approaches. If you know of any, could you please list them here?


title: ascclogit cross elasticity

text: I'm calculating elasticity right now for alternative-specific conditional logit asclogit model. Stata does not have any command or menu for elasticity calculation. I used the popular formula: (1-p)(beta)(X) for direct elasticity. But for cross elasticity is different. provided formula is: -p(i)p(j).beta where i and j are alternatives of level of dependent variable for instance in transport car and rail.


title: Calculate number of standard deviations separating two multivariate Gaussians?

text: Given a set of multivariate Gaussian distributions (from fitting a Gaussian mixture model) I would like to be able to calculate the likelihood that a data point drawn from one Gaussian will improperly be classified as belonging to one of the other Gaussians. Specifically, I have this data (ignore the inset)![data, fit to a Gaussian mixture model](http://imgur.com/G7LeWBN.jpg). If you're interested, this is figure 5 from this paper. The ellipses are drawn 3 and 4 standard deviations around each cluster. Given a model of the cluster covariance matrix as a function of the centroid, I'd like to calculate a set of centroids such that the probability of misidentifying a data point from that cluster is less than some specified probability. Put another way, I'd like to move these clusters together until the 4 standard deviation contours touch. It seems like the Bhattacharyya distance measures what I want, or is very close to it, but I'm not sure how to interpret that distance in terms of a likelihood of a point from one distribution being misidentified as a different one.