title: Meaning of a "portmanteau test"?

text: From Wikipedia > A portmanteau test is a type of statistical hypothesis test in which the > null hypothesis is well specified, but **the alternative hypothesis is more > loosely specified**. Tests constructed in this context can have the property > of being at least moderately powerful against a wide range of departures > from the null hypothesis. ... Use of such tests avoids having to be very > specific about the particular type of departure being tested. Then there are two examples   1. > In **time series analysis** , two well-known versions of a portmanteau test are available for testing for autocorrelation in the residuals of a model: the Box–Pierce test, and the Ljung–Box test. Are other tests such as the turning point test, the difference sign test, and the rank test, also considered as portmanteau tests? In Brockwell and Davis' Introduction to Time Series and Forecasting, why are they not written as portmanteau tests, but only the Box–Pierce test, and the Ljung–Box test are?   2. > In the context of **regression analysis** , a portmanteau test has been devised .... In regression, what are some examples of tests considered as portmanteau tests and not considered as such? Thanks and regards!


title: Is there any intuitive meaning to the quantity P(A|B)P(B|C)?

text: I've been playing around a little with basic probability and I was wondering if there is any intuitive meaning to either, $$ P(A|B)P(B|C) \quad\text{or}\quad \sum_{B} P(A|B)P(B|C). $$ After playing with the definitions I didn't see any obvious simplifications, but I feel there should be one. I thought it might simplify to $P(A|C)$ but that doesn't work in general. And it's definitely not $P(A,B,C)$ or $P(A,C)$.


title: Soccer game outcome, are chances equal between win/loss/tie?

text: In soccer (or any discipline with similar scoring rules), is there statistically equal chance for win, loss and tie? So, there are 3 possible outcomes, but are all 3 outcomes equally possible? This does not take into account better/worse teams, so teams are equal for the chances... So ranking, player setup or physique, etc, should not be taken into account. Also, neither should any external condition (quality of playing field, weather, ...). The reason I thinks so, is because a win is everything where scoreA > scoreB, which seems (intuitively) more likely to happen than scoreA = scoreB. If scoreB is 1, then a tie is only possible when scoreA is also 1. However when scoreB is 1, a win for A is any possible combination where scoreA > scoreB. Which is much more than the 1 single option for a tie...


title: Why is the marginal liklihood term(the denominator) in Bayes' updation not a function of theta(parameters)?

text: ![The snippet from the book is pasted below. I don't understand how the denominator is not a function of theta as it is integrated over theta](http://i.stack.imgur.com/0J1bC.png) The snippet from the book is pasted above. I don't understand how the denominator is not a function of theta as it is integrated over theta. Thank you.


title: Application of Unsupervised methods to scarce output data

text: I want to apply the Machine learning to a problem which has 9000 input rows in 132 predictive variables but outputs only for 28 of them. Should I be better of using unsupervised learning methods?. Actually the inputs are continuous while output is taken once every day.


title: Predict Variables in Multivariate Time Series

text: Please note that I am asking for advice regarding general approaches to analysis in this question. I have a time series of some $p$ variables and $n$ observations (going back several years and with resolution of approximately 15 minutes). The variables are the numbers of product $p_i$ sold, and each sale has its own observation; not all products are sold at the same time and some products may frequently be sold together. It is as yet unclear if any products are only sold conditional on some other product also being sold. If it is helpful we can think of my data a $p$ x $n$ matrix whose rows are vectors with $p$ entries representing sales per 15 minutes. I understand that I can separately analyze each of my $p$ products using a variety of methods and forecast future sales separately. However, in my (admittedly beginner's) research I have seen no discussion of how to develop a model that can predict all $p_i$ at once. Is there such a technique? Or should I form $p$ models separately?


title: What is a meaningful signifcance level when considering a pre-test and post-test model to assess learning in a college course?

text: When using a pre-test/post-test assessment to measure college learning, what is a meaningful significance level? I realize that p<.05 is the minimum criterion, yet with a correlated-groups t-test, it doesn't take much improvement to achieve this significance level. I'm getting p<.0001 by achieving only a moderate overall increase in scores with N=35.


title: What is a good measure (or set of measures) for the difference between two sample sets?

text: I'm doing some software testing where we are measuring specific latencies. Generally we run the same test several times to just eyeball the results and make sure that they are consistent across runs. I usually do this by plotting the cumulative distributions for each run together on a graph, and look for an anomalies. What more quantitative methods exist for doing this sort of thing? I would ultimately like to be able to run my tests multiple times, and determine some sort of error bound, so I can say something like "I change variable X and the median time changed Y, but that is not significant". Apologies for any vague handwaving above, I am not a statistician :)


title: Personality traits as fixed vs. random effects in a repeated measures mixed model

text: We performed a randomized controlled trial. I am interested in anxiety levels (dependent), treatment (meds vs. ultra-dance therapy), and personality variables (agreeableness, neuroticism, weeaboo). We measured anxiety levels four times. I am interested in exploring whether participants in one treatment experienced greater reductions in anxiety, and whether the personality variables are associated with reductions in anxiety over time. My understanding is that subjects are treated as random effects through the repeated measures aspect of the analysis. It seems fairly clear to me that treatment type is a fixed effect. Are the personality variables (scored 0-50 each lets say) best conceptualized as fixed, or random effects? Here is the syntax I have from SPSS if I put the personality factors as fixed effects:               MIXED ANXIETY BY TREATMENT TIME WITH PERSONALITYXXXX PERSONALITYYYY       /CRITERIA=CIN(95) MXITER(1000) MXSTEP(100) SCORING(1) SINGULAR(0.000000000001) HCONVERGE(0,          ABSOLUTE) LCONVERGE(0, ABSOLUTE) PCONVERGE(0.000001, ABSOLUTE)       /FIXED=TREATMENT TIME TREATMENT*TIME PERSONALITYXXXX0 PERSONALITYYYY| SSTYPE(3)       /METHOD=REML       /PRINT=DESCRIPTIVES  SOLUTION TESTCOV       /RANDOM=INTERCEPT | SUBJECT(id) COVTYPE(UN)       /REPEATED=TIME | SUBJECT(id) COVTYPE(DIAG).      Here is the syntax if I put personality factors as random effects               MIXED ANXIETY BY TREATMENT TIME WITH PERSONALITYXCXXX PERSONALITYYYY       /CRITERIA=CIN(95) MXITER(1000) MXSTEP(100) SCORING(1) SINGULAR(0.000000000001) HCONVERGE(0,          ABSOLUTE) LCONVERGE(0, ABSOLUTE) PCONVERGE(0.000001, ABSOLUTE)       /FIXED=TREATMENT TIME TREATMENT*TIME | SSTYPE(3)       /METHOD=REML       /PRINT=DESCRIPTIVES  SOLUTION TESTCOV       /RANDOM=INTERCEPT PERSONALITYXCXXX PERSONALITYYYY | COVTYPE(UN)       /REPEATED=TIME | SUBJECT(id) COVTYPE(DIAG).      In SPSS I put "ID" as a subject variable and "TIME" as a repeated variable in the first window. Is this necessary/a good idea? I am concerned these variables may be treated as random more than once, if that is possible. One issue with having personality variables as a random effect is the output gets quite wacky and I get an error "The final Hessian matrix is not positive definite although all convergence criteria are satisfied." Thoughts?


title: How to model the range of variables?

text: Consider a case where instead of observing actual variable values, we are observing variable range. The possible range-values could be enumerated (i.e. could be a fixed set of values, for e.g. fixed age-group ranges) or could be observation specific (and might overlap). ~~In the former case, the variable could be treated as an ordinal variable~~. (Based on @whuber's comments, treating this variable as an ordinal variable seems an unnecessary restriction). Is anything similar possible for the latter case? Are there any other methods to analyze such data? * * * I am interested in general treatment for analyzing interval data variables - as an output variable as well as an input variable to a model. Also, if there are any exploratory data analysis methods for such variables.


title: Regression with an unknown dependent variable - estimating "likelihood" to do something

text: This probably seems like a really strange question, but let me try to explain what I want to do; hopefully it will make sense. I have a data set with a couple dozen variables, such as age, level of education, self-assessed (via a Likert scale) measures of technical ability, experience, willingness to share personal information, etc etc. Using this data, I want to estimate "likelihood of participation" in using a not-yet- developed social network. Since the site doesn't exist yet, I don't have a simple dependent variable with a "1" or "0" indicating whether they've used it or not, so I can't take the usual approach with regression and estimate what impact each variable had on their participation. Is there some way to "create" a dependent variable that models "likely vs. not likely" (or something like that), and do the regression against that? I'm really lost...would be open to any ideas. Thanks, Matt


title: Integrating Prior estimates in Simrank Model

text: I am reading SimRank paper by Jeh and Widom which tries to find the similarity between objects based on the relationships between them. Effectively, SimRank is a measure that says " **two objects are similar if they are related to similar objects**." Let us denote the similarity between objects a and b by $s(a, b) \in [0, 1]$. Following the earlier motivation, a recursive equation is written for $s(a, b)$. If $a = b$ then $s(a, b)$ is defined to be 1. Otherwise, $s(a, b) = \frac{C}{\left|I(a)\right| \left|I(b)\right|} \sum_{i=1}^{\left|I(a)\right|}\sum_{j=1}^{\left|I(b)\right|} s(I_i(a), I_j(b))$ where $C$ is a constant between 0 and 1. I was wondering about whether we can incorporate any prior information about objects(Taking a semi-supervised approach). By prior information, I mean that we have similarities between some pair of objects(which may be rough/inconsistent). I tried to change the initialization in the process but proved later that the system of equations is independent of the initialization. Any thoughts on the same. Thanks


title: Coordinate descent for binomial decision in elastic net for logistic regression

text: I'm looking to implement a coordinate descent solution to the Elastic Net (Friedman et al, 'Regularization paths for generalized linear models via coordinate descent') for logistic regression. The paper explicitly gives the coordinate update for the linear case, but I'm having difficulty understanding the solution for binomial Y. Is there any literature/has anyone implemented coordinate descent for the logistic case? That is, outside their (unparseable to me) Fortran code? Thanks


title: Usage of the Friedman test

text: I have done an experiment in which I have measured the energy required to chew three different food samples on 10 subjects with 3 repetitions per sample. I am interested therefore in finding out whether there is a difference between samples according to the energy consumed. I have performed a Friedman test by ranks followed by multiple comparison procedure (Tukey's honestly significat difference for ranks) as Friedman test was significant. However someone told me that for using that test I should make sure that the error for order effect was not significant, ans I don't fully understand this. I provided my samples in the same order for all subjects. That is A,B,C -A1 B1 C1 and A3 B3 C3 for each subject. I will be very grateful if you can tell me   1. whether the test I used (Friedman) is the correct for my purpose?   2. whether the sequence in which I gave the samples affects my result from Friedman test?   3. if what I have done is incorrect, what to do instead?


title: Which k-fold cross-validation strategy is better?

text: I have met two k-fold cross-validation strategy. Both method divides the data set to k-folds, say k = 10. First one chooses the validation set elements randomly while the second one chooses incrementally. Let me explain this with an example: Method 1 chooses 3 random folds in order to use as validation set and remaining 7 folds are used as training set. And repeats this operation round count times, which is a parameter whose value is smaller than C(10, 3).               Method 1: RoundCount x ( (3 x random folds) + (7 x remaining folds) )      Method 2 takes round count as 10 so that loops for ten times. In each step, takes a different fold as validation set and the rest are used as training set.               Method 2: 10 x ( (1 x incremental fold) + ( 9 x remaining folds) )      The second method guarantees that every single sample have existed in both validation and training set. Is there an obvious advantage one of these methods to the other one?


title: Shares on total as dependent variable

text: I would like to ask for your help with the following issue. I am trying to estimate the model where dependent variable is a share of total (e.g. share of the US economy in terms of GDP on total world GDP). I use the panel data for 27 cross sections (countries) and 8 time dimensions (years). By definition, for each year the values of dependent variables across countries sum up to one. The distribution of the dependent variable is highly positively skewed with lots of observations being close to zero (or even zeros). Usually, authors propose logistic transformation of the dependent variable in terms of $y*=\ln(y/(1-y))$ for y being the original dependent variable. However, even after such a transformation the dependent variable remains highly positively skewed. The pooled OLS regression or random-effect panel regression delivers residuals that are not normally distributed and again highly positively skewed. The possible option, I guess, is to go for beta distribution of dependent variable, but not sure whether appropriate for panel data and what type of model then to use to estimate the regression (would rather have random-effect panel regression as by nature of the independent variables the fixed-effect model is not the preferred one). Thank you for any suggestions how to proceed in this case.


title: Observed versus Synthetic data

text: I am looking for studies that compare different spatial interpolation methods for observed data. However I am looking for studies that have also compared observed with generated synthetic data. For example, the comparison of spatial interpolation methods using both the data obtained from a couple of stations for one day, than the comparison with the results obtained if we replace the observed data with some generated synthetic data (random following the distribution of the observed dataset).


title: Implement Forward, Backward, Step and LASSO in VB .NET

text: My client wants me to implement Variable selection methods i.e. Forward, Backward, Step and LASSO in VB .Net platform including p-value and AIC. I have no idea about the steps involved to calculate above variable selection method. Can anyone tell me each steps including sub-steps calculation or formula. Also I am unaware about significance of P-value and AIC.


title: Differences in edition of advanced statistics books

text: In general I'm wondering if there's any danger in studying from older edition statistics books. The undergrad material seems to have remained unchanged throughout the last century, but I'm looking for 'the next step', which makes me wonder if older editions will be modern enough. For example, Applied Multivariate Statistical Analysis, Johnson/Wichern seems to be a book very much in my taste and I wouldn't mind saving $60 on a used copy from 1992. Would this have any drawbacks when self studying?


title: How to calculate 95% CI for a random effect?

text: The R code "intervals()" gives confidence intervals for fixed effects only in a mixed model. *Is there a reason why only fixed effects' confidence intervals are provided? *Is there any way to get confidence intervals for random effects as well?


title: Motivation of Projection Pursuit

text: Given I understand why PCA aims to maximize variance, I wonder why finding non-gaussian projections is interesting for Projection Pursuit, i.e., what is the motivation of this reasoning?


title: CMH-test, multiple contingency tables, unpaired entries

text: I have frequency data for 9 populations, divided into two groups, one of five populations and one of four. I am interested in whether these two groups are different. A classic way to analyse such data, I think, is a CMH test where, in my example, I could construct a 2 x 2 x k test, where k is 4 = number of pairs (i.e 4 pairs based on four population from each group). Clearly, I'd be leaving a one population out to generate a 4 x 2 x 2 table from my nine populations. While the CMH-test seems clear and straightforward for paired data and paired groups, I don't have pairs. Instead, there are numerous combinations of 4 pairs that can be generated from the 9 populations; in fact, 120. Thus, I end up, considering all 2 x 2 x 4 tables, 120 chi-square statistics, odds-ratio's and p-values. The question is how might I summarize these outputs into "one" summary. Perhaps a more formal question is "What is the distribution of a set of Chi- square statistics"? Or, is there an alternative?


title: General linear hypothesis test statistic: equivalence of two expressions

text: Assume a general linear model $y = X \beta + \epsilon$ with observations in an $n$-vector $y$, a $(n \times p)$-design matrix $X$ of rank $p$ for $p$ parameters in a $p$-vector $\beta$. A general linear hypothesis (GLH) about $q$ of these parameters ($q < p$) can be written as $\psi = C \beta$, where $C$ is a $(q \times p)$ matrix. An example for a GLH is the one-way ANOVA hypothesis where $C \beta = 0$ under the null. The GLH-test uses a restricted model with design matrix $X_{r}$ where the $q$ parameters are set to 0, and the corresponding $q$ columns of $X$ are removed. The unrestricted model with design matrix $X_{u}$ makes no restrictions, and thus contains $q$ free parameters more - its parameters are a superset of those from the restricted model, and the columns of $X_{u}$ are a superset of those from $X_{r}$. $P_{u} = X_{u}'(X_{u}'X_{u})^{-1} X'$ is the orthogonal projection onto subspace $V_{u}$ spanned by $X_{u}$, and analogously $P_{r}$ onto $V_{r}$. Then $V_{r} \subset V_{u}$. The parameter estimates of a model are $\hat{\beta} = X^{+} y = (X'X)^{-1} X' y$, the predictions are $\hat{y} = P y$, the residuals are $(I-P)y$, the sum of squared residuals SSE is $||e||^{2} = e'e = y'(I-P)y$, and the estimate for $\psi$ is $\hat{\psi} = C \hat{\beta}$. The difference $SSE_{r} - SSE_{u}$ is $y'(P_{u}-P_{r})y$. Now the univariate $F$ test statistic for a GLH that is familiar (and understandable) to me is: $$ F = \frac{(SSE_{r} - SSE_{u}) / q}{\hat{\sigma}^{2}} = \frac{y' (P_{u} - P_{r}) y / q}{y^{t} (I - P_{u}) y / (n - p)} $$ There's an equivalent form that I don't yet understand: $$ F = \frac{(C \hat{\beta})' (C(X'X)^{-1}C')^{-1} (C \hat{\beta}) / q}{\hat{\sigma}^{2}} $$ As a start $$ \begin{array}{rcl} (C \hat{\beta})' (C(X'X)^{-1}C')^{-1} (C \hat{\beta}) &=& (C (X'X)^{-1} X' y)' (C(X'X)^{-1}C')^{-1} (C (X'X)^{-1} X' y) \\ ~ &=& y' X (X'X)^{-1} C' (C(X'X)^{-1}C')^{-1} C (X'X)^{-1} X' y \end{array} $$   * **How do I see that** $P_{u} - P_{r} = X (X'X)^{-1} C' (C(X'X)^{-1}C')^{-1} C (X'X)^{-1} X'$?   * **What is the explanation for / motivation behind the numerator of the 2nd test statistic?** - I can see that $C(X'X)^{-1}C'$ is $V(C \hat{\beta}) / \sigma^{2} = (\sigma^{2} C(X'X)^{-1}C') / \sigma^{2}$, but I can't put these pieces together.


title: Density function of the max of smallest and largest observation

text: Consider $n$ independent uniform random variables $X_i \sim U(-\theta,\theta)$, and let $Y_1 = \min(X_1, \ldots, X_n)$ and $Y_n = \max(X_1, \ldots, X_n)$ . What is distribution of $Z = \max (-Y_1,Y_n)$, given that the joint PDF of $Y_1$ and $Y_n$ is $$ f(y_1, y_n) = \frac{n(n-1)}{(2θ)^n} (y_n-y_1)^{n-2} , \quad -\theta < y_1 < y_n < \theta ? $$


title: ROC curve and confusion matrix in classifier performance evaluation

text: I applied two different classifiers against the same validation set. It turns out that classifier A is better than classifier B in terms of ROC curve. However, classifier B is better than classifier A in terms of confusion matrix. How to explain this kind of contradiction?


title: How to find conjugate prior for a given distribution?

text: Given a distribution $\{p(X|\theta),\theta\in\Theta\}$:   * how do you show the existence/non-existence of its conjugate prior?   * what are some general ways/principles to construct/find its conjugate prior, if there is one?    * If it exists, how do you show it is unique or not?    * If not unique, which one is the most "canonical/classical" (in your definition) one among all?


title: Determining the probability of $X_2 \ge X_1$ given they have different probability functions

text: Suppose that I have a random variable $X_1$ which is normally distributed, and a random variable $X_2$ having the density function shown in the figure below. How would I determine ${\rm P}(X_1 \le X_2)$? ![enter image description here](http://i.stack.imgur.com/ci1mv.jpg)


title: Get canonical loadings (i.e. scores) from weights in r package PMA

text: I want to perform regularized canonical correlation between two matrices with more variables than observations (same subjects), one of which is very large (~18000 columns). The only r package that could handle the matrix dimension was PMA (tried mixOmics, RGCCA...). The problem is that the output of the function CCA gives only the canonical variates and correlations, and the canonical weights. For a good interpretation of the results the canonical scores are necessary, but this is (weirdly) not provided from the CCA function. I am not very proficient with matrix algebra (and statistics in general) and I am not sure about how to compute this. Example from the package documentation:               ### create matrices     > u <- matrix(c(rep(1,25),rep(0,75)),ncol=1)     > v1 <- matrix(c(rep(1,50),rep(0,450)),ncol=1)     > v2 <- matrix(c(rep(0,50),rep(1,50),rep(0,900)),ncol=1)     > x <- u%*%t(v1) + matrix(rnorm(100*500),ncol=500)     > z <- u%*%t(v2) + matrix(rnorm(100*1000),ncol=1000)     ### perform canonical correlation (3 canonical variates)         > out <- CCA(x,z,typex="standard",typez="standard",K=3)          > print(out,verbose=TRUE)     Call: CCA(x = x, z = z, typex = "standard", typez = "standard", K = 3)          Num non-zeros u's:  59 88 75      Num non-zeros v's:  180 154 164      Type of x:  standard      Type of z:  standard      Penalty for x: L1 bound is  0.3      Penalty for z: L1 bound is  0.3      Cor(Xu,Zv):  0.9578624 0.93371 0.9418701          Component  1 :          Row Feature Name Row Feature Weight     1                 1              0.112     2                 2              0.080     3                 3              0.124        4                 4              0.165     5                 5              0.087     ........     ........                 Column Feature Name Column Feature Weight     1                    10                 0.006     2                    15                -0.027     3                    25                -0.025     4                    28                 0.030     5                    35                 0.035     ........     ........      similar for Components 2 and 3. The other doubt is about the significance of the results: are the non zero variables given already computed as "significantly different than 0"? Although a significance threshold is not asked. Thank you very much for any suggestion!


title: The probability of detecting a defect given defect growth when the two variables are independent

text: I am dealing with a real world problem and just wanted to check if I am going about this the right way. I have the probabilities that   * a defect will grow within a particular structure, say ${\rm P}(X_1)$   * a person will be able to detect a defect, say ${\rm P}(X_2)$. These are two completely independent variables. I am trying to determine the probability that a defect will be detected given that there is defect growth in the structure. Can I just multiply the two?


title: I have a problem with bayes and train lines

text: I have this problem which I can't fully understand it: > Assume you have entered in a foreign city of unknown size. At entrance you > see a tramcar with number $100$. Let us assume that the tramcar lines of the > city are numbered with integers starting from one, and that every line has > an equal probability of being spotted. What can you learn about the number > of tramcar lines on the basis of your observation. Notice that this problem > is a discrete version of the previous one. Assuming that the number of > tramcar lines $\theta$ is large, you can approximate sums with integrals and > use the following general result: > > $$\int_k^{\infty} t^{-m}\,\mathrm{d}t=\frac{k^{(1-m)}}{m - 1}$$ > > For instance, you can use a prior of the form $\theta^{-1}$ and calculate > the posterior median. Now I know how to do it with considering one line, so I will get the median if there only exist one, but multiple I'm lost, and I think is with multinomial distr. But I'm not sure how to set it up, I have read the German tank problem which it sounded similar and the explanation was very intuitive and clear, I was hoping for one more or less like so, I'm not an statistician, more a programmer dealing with bayes atm. Any help? Thanks.


title: What is the difference between accuracy and agreement?

text: According to Manning et al. (p. 155) **accuracy** is the sum of the diagonal in the confusion matrix divided by the sum of all items. On the other hand, following Artstein and Poesio (p . 558) precisely this seems to be also called **agreement** (without chance correction). Is there any conceptual difference?


title: Trend test on historical stock prices, where to start?

text: I am a beginner on statistic but come across an urgent situation that I have to work with statistics of stock prices data. I have learned that ones can do a hypothesis test that rejects the Random Walk hypothesis in order to show that the data have trends. I have tried to understand it but that seems impossible for me as a completely beginner. (beside simple concepts like probability and linear regression, I have just learned about the hypothesis test along the way!) I have tried finding if there are easier alternative ways that I can understand but am feeling lost as there are too many terms that I don't know such as parametric-nonparametric thing, t-statistics, residuals etc.. Not to mention that I can't randomly pick a test strategy to do the work, could I? In fact, I am not doing any deep analysis at all. All I want to do is to prove that a subset of my data is statistically increasing or not (in fact it is, visually), as follows. **I have a time series of historical prices, $\{p_0, p_1, p_2,...,p_n\}$. And I have an algorithm to pick 2 points $p_{t1}$ and $p_{t2}$ (with $ {t2} > {t1}$) such that $p_{t1} < p_{t2}$. Then I want to show that the prices between $p_{t1}$ and $p_{t2}$ is statistically increasing.** Please guide me to do that or advise a learning path that leads straight to the topic. Thank you. **Additional Information** The picture below shows examples of the case where the prices between the two points are increasing and the case where not. ![Example of the case where the prices between the two points are increasing and the case where not.](http://i.stack.imgur.com/MNkQH.png) My algorithm involves the use of MACD and fuzzy logic but it's just a trial- and-error thing and may looks basically like nothing more than asking a person to pick two points as mentioned. My intention is that the algorithm chooses $p_{t2}>p_{t1}$ due to that the data have the increasing trend from $t1$ to $t2$, not by chance


title: Mixing continuous and binary data with linear SVM?

text: So I've been playing around with SVMs and I wonder if this is a good thing to do: I have a set of continuous features (0 to 1) and a set of categorical features that I converted to dummy variables. In this particular case, I encode the date of the measurement in a dummy variable: There are 3 periods that I have data from and I reserved 3 feature numbers for them: 20: 21: 22: So depending on which period the data comes from, different features will get 1 assigned; the others will get 0. Will the SVM work properly with this or this is a bad thing to do? I use SVMLight and a linear kernel.


title: Any thoughts in PLS as an alternative to SEM if you have a small sample size?

text: I have put a question on here before re SEM and small sample sizes. I have 3 factors each with 3 manifest variables each. I only have 90 records. I want to find out the amount of covariance between the factors AND between ALL the manifest variables. I also want to determine whether the 3 factors predict a forth performance measure. I was wondering whether Partial Least Squares (PLS) would be an alternative? Does this make sense? Any thoughts?


title: Confusion about the EM algorithm

text: I am reading through the EM-Algorithm but on the slide 39, I don't get how $$P(D) = P(D|A)P(A) + P(D|B)P(B)$$ I am trying to understand it in order to get my head around modelling with mixture of Gaussians. Can anyone explain me how to get that result?


title: What if only control variables are significant in a differences-in-differences analysis?

text: Regarding the standard DID model: $$ y=\alpha+\beta_1\text{treat}+\beta_2\text{post}+\beta_3\text{treat⋅post}+u $$ What exactly does it mean if say $\beta_3$ is not statistically significant, but $\beta_1$ is? Does the significance of $\beta_1$ just mean that my control group and my experimental group differed in the very beginning in a statistically significant manner? In addition, say that you add in additional control variables. Does it matter if the coefficients on any of those control variables is significant, or do we only care about the control variables in how they affect the value and significance of $\beta_3$? Essentially: how do you interpret the coefficients on your control variables if they end up being significant in a DID regression?


title: Data set to experiment with multivariate count regression

text: I have a model for predicting multiple count variables (multivariate count regression) given some covariates. Are there any publicly available datasets I could experiment with?


title: Consequences of high power for our belief in the null hypothesis?

text: Suppose that we choose the size of the test to be $\alpha = 0.05$, and based on our sample size $n$ and magnitude of the absolute value of the estimate, we determine that the test's power (i.e., 1 - Prob(Type 2 Error)) is $99.9\%$. If I then fail to reject $H_0=0$ against $H_a \not= 0$, what consequences does this have for our belief that $H_0$ is correct, assuming that we have the correct null hypothesis distribution of the test statistic? What if the power was $50\%$?


title: How to calculate median of distributed data?

text: How to find the median (or approximate median) of n chunks of values. The values in each chunk can be assumed to have been sorted if it helps. In addition, the existence of a histogram for each chunk is allowed.


title: What is the distribution that can properly describe the PE fluctuation of a stock

text: I have observed the historical PE (price / profit) value of a stock and realized that it roughly follows a log normal distribution. However, even when the next earning data point is easily predictable, this distribution cannot be used to predict the distribution of next stock price data point, because, the next PE value depends on the previous one, thus it is not set of repeated independent trials. In other words, for a regular random variable that follows the log normal distribution, the next value can be anything, regardless of the current value. But for the PE value of a stock, it must change gradually. It is also not exactly a random walk, because for the random walk, the next step's direction is also independent of the previous step. This is not true for the PE, when it is too low, it is much more likely to go up than go down. So my question is, what kind of random process can properly model the movement of the PE, which changes gradually, but overall, follows a log normal distribution?


title: valid instrument for oil consumption in IV model

text: I want to run a gdp vs. oil consumption model where oil consumption is suspected to be endogenous - correlated with the error terms. Can a variable correlated with world oil price but not with the gdp of the country of interest be a valid instrument?


title: Comparing regression coefficients across groups in multinomial logistic regression

text: Suppose that I fit a Bayesian multinomial logistic regression model where the dependent categorial variable indexes $x$ groups, and the predictors are the same across groups. I now have $x - 1$ sets of $p$ regression coefficients and want to compare the regression coefficients of one of the groups to the two others. One simple ad hoc way to do this is to compare the posterior distributions of the coefficients. But I'm sure there is a better way. Thanks.


title: When interpreting interactions in a factorial ANOVA, is it necessary to look at the residual cell means?

text: In their 1989 paper "Definition and Interpretation of Interaction Effects", Rosnow & Rosenthal write: > When interaction is claimed in a factorial arrangement, the results almost > always require more detailed analysis than is typically reported in our > primary journals. In reporting interactions, research psychologists have > gotten into the habit of examining only the differences between the original > cell means (the simple effects) instead of more properly examining the > residuals, or leftover effects, after the lower order effects have been > removed. And then later in the same article: > The point of this article is to emphasize that if investigators are claiming > to speak of an interaction, the exercise of looking at the corrected cell > (or condition) means is absolutely necessary. Is this really true? The reason I ask this question is because it seems to me that despite this article being published 25 years ago, it seems to me that people generally still generally interpret interactions in the way condemned by Rosnow and Rosenthal. Is that because a bad habit has persisted? Or is it the case that Rosnow and Rosenthal were wrong in the first place?


title: Sufficient statistic with 4 dimensions

text: I have an i.i.d. sample with this density: $ f_\theta (x) = k * exp (- (x - \theta)^4 ) $ How can I obtain the sufficient statistic? Is the statistic that I find also complete? Thank you very much!


title: Is it possible to apply KL divergence between discrete and continuous distribution?

text: I am not a mathematician. I have searched the internet about KL Divergence. What I learned is the the KL divergence measures the information lost when we approximate distribution of a model with respect to the input distribution. I have seen these between any two continuous or discrete distributions. Can we do it between continuous and discrete or vice versa?


title: Fitting a constrained natural spline

text: **Why I think I need to fit a constrained natural spline** I'm working with IBI (interbeat intervals). The measured IBI's include some invalid values. I'm working on an imputation algorithm to impute for the invalid IBI's. Suppose I have an interval with 5 invalid IBI's. Part of my idea on imputing for the invalid IBI's is to have a natural spline fit that is constrained between the last valid point before these invalid IBI values and the first valid point after these invalid points. **Is there any way that I can fit a natural spline constrained to pass through a certain pre-specified points?** I am working in R, and R implementations would be particularly useful. As an example, consider the data below:               data <- data.frame(X = rep(NA, 1000), Y = rep(NA, 1000))     set.seed(123456)     data$X <- cumsum(abs(c(0, rnorm(999, 0.8, 0.05))))          data$Y <- c(rnorm(100,700, 50), rnorm(100,800, 50), rnorm(100,900, 50), rnorm(100,800,         50), rnorm(100,700, 50), rnorm(100,800, 50), rnorm(300,700, 50), rnorm(100,600, 50))      Now, suppose that I'd like to have a natural spline with knots on every 20 seconds (assume time is in second) AND to have the fitted natural spline fit to pass through obs 4 and obs 15 (row # 4 and row # 15). Is there anyway to do that? **My solution** : one solution is to use very close knots in those points.


title: Logistic regression for abiotic influences on behavior

text: What I am looking to do is test for a correlation between an activity (in this case nesting) with cumulative rainfall from the previous two weeks. For example, say one individual nested on DayX where the total rainfall from the previous 2 weeks is 4cm and another individual nested on DayY with 8cm of prior rain, and there are days where it had rained with no nesting events. An example of the data is such:               Day 2wk_rain number_nested number_available     1     2            0              9      2     2            0              9     3     4            0              9     4     3            1              9     5     6            1              9     6     2            2              9     :     :            :              :     15    8            9              9      Analyzing this with cumulative rainfall throughout the period seemed straightforward, but when I decided to look further at a more localized temporal scale it threw me for a loop. My main objective is to find out whether rainfall can trigger nesting behavior within a few days, hence scaling it down. I still think I need to stick with logistic regression due to the nature of the response variable, but how I go about analyzing it and presenting it in an understandable way seems to elude me.


title: Conceptual questions on Estimation theory: Beginner level

text: I have just started to study Estimation theory and following the book "Fundamental of Statistical Signal Processing" by Steven M. Kay Vol 1. Most of the examples are based on the assumption that data follows a Gaussian distribution. COnsider, the example of estimating an unknown parameter $A$ from the relation: $x[n] = A + w[n]$ where w[n] is said to be $i.i.d$ random gaussian noise of known variance and zero mean. I have few conceptual questions and would appreciate intuitive notion behind them Q1: Performing, $x[n]-w[n] = A$ would yield the unknown parameter A since w[n] is known. Then why does one need to consider the estimation using the estimator $\hat{A} = \frac{1}{N}\sum_{n=0}^{N-1}x[n]$ (Pg 12, Sec 2.3 Unbiased Estimators) Q2: In general, does the number of samples, n play any role in estimation? Say n= 10 and n=1000, would this consideration affect the possibility of finding better estimates? Q3: In general, why do we consider the samples and the noise to be $i.i.d$ ?? Is it ever possible that the samples are not iid? If so then when are they correlated? Thank you


title: Estimation of a state-space model using Bayesian analysis with the Metropolis-Hastings algorithm

text: I have the following state-space model: $$\begin{aligned} y_t&=c+Ax_t+q_t, &q_t \sim \mathcal N(0,Q), \\ x_t&=\mu+Bx_{t-1} + v_t, &v_t \sim \mathcal N(0,R), \end{aligned} $$ where the series $x_{t=1}^T$ is latent. I am interested in estimating $\Phi=[ x_{t=1}^T, A,B,R,Q,c,\mu ]$ via Bayesian analysis using the Metropolis-Hastings algorithm. That is, I want to estimate the following distribution: $$ f(\Phi|y^t) \propto f(y_t|\Phi)f(\Phi). $$ $f(y_t|\Phi)$ is a maximum likelihood function which can be constructed with the Kalman Filter algorithm. I am unsure about the second function: $f(\Phi)$ (the prior). What is the proper prior for the whole series $x_t^T$? Also, as $f(\Phi)=f(x_t,A,B,R,Q,c,\mu)$ I am not sure if it makes sense to have independent priors as $\mu$ and $R$ are not independent of $x_t^T$.


title: How to interpret coefficients in the spatial lag and spatial durbin model?

text: How do I interpret the coefficient for the spatially lagged dependent variable in the spatial lag model or spatial autoregressive model (SAR) (y = αιn + ρWy + Xβ + ε) and how do I interpret the coefficient of the spatially lagged dependent variable and the coefficients of the spatially lagged independent variables in the spatial Durbin model (SDM) (y = αιn + ρWy + Xβ +WXθ + ε)? I am aware of the interpretation of the direct and indirect effects for the independent variables, but what is with the other coefficients. Most articles I read using the SDM just did not interpret the coefficient of the spatially lagged dependent variable and the coefficients of the spatially lagged independent variables at all. Thanks for your help!


It seems that KNN is a discriminative learning algorithm but I can't seem to find any online sources confirming this. Is KNN a discriminative learning algorithm?


title: Comparing polynomials

text: I've got a bunch of data on pop singers' performance on the Hot 100 charts over time, and I'm trying to compare the early part of different artists' careers. For example, I might look compare Miley Cyrus' career to the commensurate period of Cher's career: ![enter image description here](http://i.stack.imgur.com/3XWQb.png) (The y axis is total "top 10 points" -- 10 points for a weekly #1, etc.) I'm trying to identify the singer among the 2000-some in the database who had the most similar early career. Running a simple correlation on the first 7 years of everyone's career wasn't very fruitful: ![enter image description here](http://i.stack.imgur.com/iHDJB.png) My next idea was to fit a binomial to each person and compare. Here's Beyonce vs. The Beatles through six years, for example: ![enter image description here](http://i.stack.imgur.com/992us.png) The trouble is, I have no idea how to compare a bunch of polynomials to figure out where they're similar. Any ideas how to do that, or otherwise how to compare these datasets?


title: Offset in Poisson regression

text: Consider an offset term in a Poisson regression: $$\log \mu_x = \log t_x+ \beta_0 + \beta_{1} x$$ To interpret $\beta_0$, would you need to consider $(\beta_0+ \log t_x)$? Because what if $t_x \neq 1$? Is the interpretation of $\beta_0$ the mean number of events when $t_x = 1$ and $x=0$? Suppose $t_1 = 2, t_2= 3$ and $t_3 = 4$. Then there is no $x$ such that $t_x = 1$. Also $x \neq 0$. Also is the interpretation of $\beta_1$ the following: The mean number of events comparing $x+1$ and $x$ for a fixed $t_x$?


title: Conditional fixed effects Poisson *with population exposure*

text: Reading a paper, the authors use a "Conditional fixed effects Poisson with population exposure." I wasn't sure what "population exposure" meant, so I looked it up. It seems like the term comes from epidemiology (not my field). Could I check to see if I've understood it correctly? Poisson models try to estimate a rate of event occurrence within a population. The relevant population is subject to interpretation. In epidemiology, the population might be called the "exposure" group -- e.g. those exposed to the disease risk. So, by saying "population exposure" a Poisson model with population exposure calculates the event rate as: $\ \ \ \ \ \log(\frac{\text{event_ count}}{\text{region_population}}) = X\beta$ $\Rightarrow \log(\text{event_ count}) - \log(\text{region_population}) = X\beta$ $\Rightarrow \log(\text{event_ count}) = X\beta + \log(\text{region_population})$ I appreciate any help you can give toward understanding this.


title: Given the distribution of a variable X, what is the distribution of f(X)?

text: Suppose I know the mean and standard deviation of a (roughly) normally distributed variable x. Is there any way for me to calculate the mean and standard deviation of f(x)? The particular problem I am working with involves needing to find the mean and standard deviation of sine and cosine of variables whose mean and standard deviation I already know.


title: variational bayes vs expectation propagation

text: Just wondering what are the advantages of using one over the other. I'm just looking for some general answers here. For starters: * VB gives a guaranteed lower bound for the likelihood.  * EP is faster? VB has quite a lot of iterations till convergence, unless someone came up with a faster method. * With EP you have to specify the family of the approximate posterior


title: Include dummy variables in model - always necessary?

text: I have some survey data collected from 15 government organizations. Using this data, I would like to run an OLS regression to estimate the influence of variables x1 and x2 on y. I have run two models: one including 14 dummy variables accounting for potential differences from organizational membership, and one not including any dummy variables. The difference between R-square values between the models is small (less than 0.02), and the adjusted R-square difference much smaller. Moreover, the coefficients for x1 and x2 and significance levels are not substantially different between models, and none of the dummies has a statistically significant relationship with the dependent variable. I am hoping for some insight into the following questions:   1. From a theoretical standpoint, is it necessary to include dummy variables for organizational membership (in order to protect against omitted variable bias or some other reason)? In other words, must I include the dummy variables on principle?   2. If I must include them on principle, is there any post analysis test that I can use to claim that the model without the organizational dummies should be used in the interpretation of the results? (Something similar to the Hausman test for distinguishing between fixed and random effects models, for instance.) Please let me know if I can provide any more information that will help.


title: How to decompose the probability of sum of a series random variables？

text: Let $A_1, A_2, \ldots$ be independent random variables with known distributions. I'm wondering if it's possible to equivalently decompose the probability above into separated constraints?   For example:   \begin{align*} \Pr(A_1 < b_1) & > 1-p_1 \\ \Pr(A_2 < b_2) & > 1-p_2 \\ &... \end{align*}


title: Image vectorization in matlab

text: I need to convert images into vectors (image vectorization). I have 165 images in total, divided into 15 subjects of 11 images each. Following is code i have written to convert images into vectors.                    no_of_images = 165;                   class=dir('<path of the images>*.pgm');         for i=1:length(class)             file_name=strcat('<path of images>',class(i).name);             a = imread(file_name);         end                [m, n] = size (a);                        data = (zeros(m*n, no_of_images));                      % Create data matrix     for i = 1 : no_of_images         class=dir('<path of the images>*.pgm');         for j=1:length(class)             file_name=strcat('<path of the images>',class(j).name);             image = imread(file_name);         end              data(:, i) = reshape (image, m*n, 1);     end;      When i run this code, i goes into an infinite loop of only zeros. Where am i going wrong?


title: Generating a correlated data matrix where both observations and variables are correlated

text: I am trying to generate a simulated data matrix that is correlated by both observation and variable directions. So far I know how to do this for **variable x variable**. ![enter image description here](http://i.stack.imgur.com/NgP8j.jpg)                # correlated matrix between variables          n = 200         p = 100          CRMt <- matrix(NA, nrow = p, ncol = p)         diag(CRMt) <- 1         CRMt[upper.tri (CRMt, diag = FALSE)] <- 0.5         CRMt[lower.tri (CRMt, diag = FALSE)] <- 0.5          L = chol(CRMt)# Cholesky decomposition     p = dim(L)[1]          set.seed(999)     M = t(L) %*% matrix(rnorm(p*n), nrow=p, ncol=n)     M1 <- t(M)     rownames(M1) <- paste("S", 1:200, sep = "")      colnames(M1) <- paste("M", 1:100, sep = "")     cor(M1)      Now say I want to create a data matrix that also follows the following **observation x observation** correlation matrix. ![enter image description here](http://i.stack.imgur.com/XkDuZ.jpg)               OCRMt <- matrix(NA, nrow = n, ncol = n)     diag(OCRMt) <- 1     OCRMt[upper.tri (OCRMt, diag = FALSE)] <- 0.3     OCRMt[lower.tri (OCRMt, diag = FALSE)] <- 0.3      How can I do this ?


title: GARCH modelling

text: Is this simplistic way of describing GARCH processes correct: > Future prices do not depend on previous prices per se, but rather on the > previous variance. If this is incorrect, is there a model that does this?


title: Estimation of constant in MA(1) process using least square estimate

text: I would like to estimate the constant in MA(1) process using LSE, but am finding it very difficult.


title: unbalanced samples random Forests

text: I am trying to predict species presence or absence using randomForest in R (classification). In fact, I am trying to do it for several species, in separate models. For a couple of the species, the training data are quite unbalanced e.g., 70 observations of species presence, and and 6500 observations of species absence. This is my code:               #read in data frame containing observations of species presence/absence and predictor     variables     mydata <- read.csv('mydata.csv')          #fit random forests model     fitmodelA <- randomForest(SPECIESA ~ var1 + var2 + var3 + var4 + var5 +var6 + var7 +   var8 + var9 + var10, data=mydata, mytry=3, ntrees=500, replace=TRUE, importance=TRUE,   keep.forest=TRUE)          #predict to new data     predictmodel <- predict(fitmodelA, newdata, type="prob")      In the output prediction, almost the entire study area is predicted with prob > 0.7. I take this to be predictions of species occurrence? or is it the probability of species absence? I want to try to balance the data by forcing the model to select equal sample sizes from observations of presence and absence, e.g., adding the argument               sampsize(70,70)      But I get the error message "Error in if (ncol(x) != ncol(xtest)) stop("x and xtest must have same number of columns")" What am I doing wrong here?


title: Continuous AR model from a Discrete AR model

text: I have a model of an $\mathrm{AR}\left(2\right)$ process, thus: $X_n + a_1X_{n-1}+a_2X_{n-2}=Y_n \qquad \mbox{for} \qquad n=0,\pm1,\pm2,\ldots $ and I have the equivalent stochastic differential equation in Continuous Time: $\ddot{X}\left(t\right)+\alpha_1 \dot{X}\left(t\right) + \alpha_2 X\left(t\right) = \varepsilon \left(t\right)$ But I now want the $\mathrm{AR}\left(2\right)$ process to have a seasonal term aswell. How do I write the $\mathrm{AR}\left(2\right)$ process, and the above stochastic differential equation, to include this?


title: estimate a distribution parameters only by data mean and std. dev

text: I need to estimate a truncated gamma distribution parameters (shape , scale). But, I only know the data mean and std. dev. I do not know the data set. Given the mean and std. dev. of a data set from a truncated gamma distribution, how to find shape and scale of the distribution parameters ? I know MLE may be useful for solving this problem. But, they depend on knowing the whole data set. Any help would be appreciated. * * * **UPDATE** I have used the Method of moments here http://en.wikipedia.org/wiki/Method_of_moments_%28statistics%29#Example But, the QQ-plot from the estimated parameters (got by the above method) is very bad. I think this is because the method of moments do not consider the truncation points. * * *


title: Angled brackets in notation (expectation or inner product doesn't seem to make sense)

text: I am reading an internal paper at my university that refers to volatility (for a variable not explicitly mentioned) as: $\langle X,X\rangle = \int_0^T \sigma_t \, \mathrm{d}t.$ I am unfamiliar with expectations using this notation with commas, and subsequent uses do not indicate inner products, inner products wouldn't really make sense, and the $\langle \rangle$ are used rather than the $<>$ as well. Volatility of volatility is referred to as: $\langle \sigma^2, \sigma^2\rangle_T$. Thanks for any assistance!


title: Statistical test for feature selection

text: Suppose I have a procedure to select K features out of M. I repeat the procedure N times on, say bootstrapped datasets, and count how many times feature #1 appeared among selected, denote k1, how many times feature #2 appeared, denote k2, etc. Question 1: What would be a suitable statistical test to check from the counts k1,k2,...,kM, whether the selection procedure is random or not? Question 2: What would be a suitable statistical test to check from ki/N, whether feature i appears ki/N times by chance, or appears more/less frequently than it would simply by chance? Question 3: is there a theoretically grounded way to set N, i.e. how many times to repeat the selection procedure, for the test to be conclusive?


title: Correlating continuous clinical variables and gene expression data

text: In SVM (linear kernel) classification analyses of a data-set of gene expression (~400 variables/genes) for ~25 each of cases and controls, I find that the gene expression-based classifiers have very good performance characteristics. The cases and controls do not differ significantly for a number of categorical and continuous clinical/demographic variables (as per Fisher's exact or t tests), but they do differ significantly for age. Is there a way to show that the classification analysis results are or are not influenced by age? I am thinking of reducing the gene expression data to principal components, and doing a Spearman correlation analysis of the components against age. Is this is a reasonable approach? Alternately, can I check for correlation between age and class-membership probability values obtained in the SVM analysis. Thanks.


What does it mean for a study to be over-powered? My impression is that it means that your sample sizes are so large that you have the power to detect minuscule effect sizes. These effect sizes are perhaps so small that they are more likely to result from slight biases in the sampling process than a (not necessarily direct) causal connection between the variables. Is this the correct intuition? If so, I don't see what the big deal is, as long as the results are interpreted in that light and you manually check and see whether the estimated effect size is large enough to be "meaningful" or not. Am I missing something? Is there a better recommendation as to what to do in this scenario?


title: Joint test of change in variance across matched groups

text: I am monitoring variability in sentencing outcomes (custodial sentence length) amongst matched groups. The aim is to determine whether sentence variability went up or down after a policy was implemented. I have about 40 groups of sentences, which are matched on the characteristics of the case. For example, group 1 comprises sentences for people who committed ABH, who were under the influence of drugs, and used a weapon. Group 2 comprises people who committed ABH, who were remorseful, and had no aggravating factors in the case. For each group, I can further separate sentences outcomes into two groups: those sentenced before a policy was implemented, and those sentenced after the policy was implemented. I can then compare sentenced outcomes in the 'before' and 'after' groups to see whether variability seems to have declined. I can test each group individually to see whether variability of sentence outcome went up or down, using an F test. However, the results of these tests vary, and I realise that I should really be using a joint test in any case. I'm also not getting enough power on a group-by-group basis: for most tests, we cannot reject H0, that variability did not change. This is despite the fact that variability has declined for the vast majority of groups, which suggests that if I somehow combine results from all groups, it should be possible to find a statistically significant decline in variance overall. Could anybody suggest and appropriate test or how I might proceed? I am very grateful for any help. A few clarifications:   1. I would be happy to assume that amongst a matched group (log) sentence outcomes are drawn from a normal distribution. This assumption will not be perfect, but is probably adequate.   2. The size of matched groups differs. For instance, the group 1 'before' implementation group may have 100 sentences, and the 'group 1 after implementation' group may have 75 sentences. The same figures for group 2 may be 66 and 123.


title: Finding Bayes estimator

text: Suppose $X_1,\ldots,X_n$ is a random sample of distribution $f(x|\theta)=\theta (1-\theta)^x,x=0,1,2,\ldots\ ,0<\theta<1$. How can I calculate Bayes estimator of parameter $\gamma(\theta)=\frac{1-\theta}{\theta}$ under loss function $L(\theta,\delta)=\displaystyle\frac{\theta^2}{1-\theta}(\delta-\theta)^2$ and prior distribution $U(0,1)$


title: How should I cluster my data and then run separate algotithms on each cluster?

text: I'm new to statistics and R. So, I have a dataset on which I want to build a predictive model. Since there is quite a variability in the input variables, I thought it would be better to cluster my data first (using kmeans) and then build separate models on each cluster (e.g. Random Forest). Using 10-fold cross-validation the results look very good however I'm not sure if what I'm doing is correct. Firstly, I cluster my data using kmeans but **I'm also using the dependent variable as well**. Then I fit separate models to each cluster and perform that 10-fold CV for each cluster separately. What if I had a new data point? How would I assign it to a cluster (since I won't have the value of the dependent variable)? Would it be more correct to cluster my initial data without using the dependent variable?


title: Max Entropy Solver in R

text: I am trying to solve for the parameters in the maximum entropy problem in R using the nonlinear system: $\ \int_l^u e^{a+bx+cx^2}dx=1$ $\ \int_l^u x e^{a+bx+cx^2}dx=\mu$ $\ \int_l^u x^2 e^{a+bx+cx^2}dx=\mu^2+\sigma^2$ where $l$ and $u$ represent the lower and upper bound on the support of the distribution, and $a,b,c$ are the parameters of interest. I wrote an implementation based on minimising the sum of the squared differences, and then optimising using BBoptim, but this method is very hit- or-miss, and depends heavily upon me inputting reasonable starting parameters, and quite tight lower and upper bounds on the range over which to search. I have included the code below:               MaxEnt <- function(y) {     G <- function(x) exp(y[1]+y[2]*x+y[3]*x^2) ;     H <- function(x) x*exp(y[1]+y[2]*x+y[3]*x^2) ;     I <- function(x) x^2*exp(y[1]+y[2]*x+y[3]*x^2) ;     (integrate(G, lower=l, upper=u)$value - 1)^2+(integrate(H, lower=l, upper=u)$value - mu)^2+(integrate(I, lower=l, upper=u)$value - mu^2-sd^2)^2}     BBoptim(c(-log(sqrt(2*pi)),0,0), MaxEnt, method=1)      Does anyone know of any better ways of doing what I am trying to achieve using R? Many thanks in advance. ============ EDIT JAN 17th ============ I have tried solving this problem using whuber's approach and so far this is what I get: Let $\ y= \frac{x-l}{u-l}$. Then $\ dx = (u-l)dy$. Then since $\ x= (u-l)y+l$, $\ a+bx+cx^2 = q+ry+cy^2$ where $\ q=a+bl+cl^2,r=(u-l)(b+2cl),s=c(u-l)^2$. The three original equations then become: $\ \int_0^1 e^{q+ry+sy^2}dy=\frac{1}{u-l}$ $\ \int_0^1 ((u-l)y+l) e^{q+ry+sy^2}dy=\frac{\mu}{u-l}$ $\ \int_0^1 ((u-l)y+l)^2 e^{q+ry+sy^2}dy=\frac{\mu^2+\sigma^2}{u-l}$ whuber then suggests rewriting $\ e^{q+ry+sy^2}$. Letting $\ \nu = \frac{-r}{2s}$ and $\ \tau = \sqrt{\frac{-1}{2s}}$ we can see that $\ e^{q+ry+sy^2} = e^{q-\frac{r^2}{2s}}e^{-0.5\left(\frac{y-\nu}{\tau}\right)^2}$. The problem I now have is that using whubers approach we must have that $\ \tau^2 >0 $ (?). This means I need $\ s<0$. By definition, $\ s <0$ only if $\ c <0$ also. However, I have an example where the density that has maximum entropy has $c>0$. I am probably going wrong somewhere, and would be grateful if somebody could point out where. Thanks!


title: Least squares problem

text: You are given a "blackbox": $RSS = (Y - X\beta)'(Y-X\beta)$ for any linear model $Y=X\beta+\epsilon$ You have $n=60$ observations, and two predictor variables. You want to test $H_0: \beta_1 = \beta_2$ in the model: $Y = \beta_0+\beta_1X_1+\beta_2X_2+\epsilon$ Describe how you would accomplish this using the "black box"? My professor has outlined changing $H_0: \beta_1 = \beta_2$ to be $\theta = \beta_1 - \beta_2 = 0$, so $H_0: \theta = 0$. Then rewrite the linear regression (scalar) form as: $Y = \beta_0 + (\beta_2+\theta) X_1 + \beta_2 X_2 + \epsilon$ **I don't understand where we take the algebra from here...**


title: How to visually show the convergence of a process?

text: I have some set of measurements that I have represented as vectors $x^t$ for $t \in \{ 1, 2, ...\}$. I want to test "convergence" of the process (visually) in some sense. I thought maybe I could run PCA on $x^t$ and then plot the projections to the first principal component, 2nd, etc... and see whether these plots converge to something. What is the right way to do it? There is variance in between samples, of course, so just plotting the first PCA component will not necessarily converge to something. But what should I look for, and how can I make it visually compelling? (Or maybe there is another, better way of testing the "convergence" of these $x^t$?) (Please note: the $x^t$ are not coming from an MCMC process, but maybe tools from there could be used?)


title: Questions concerning copulas

text: I am new to the topic of copulas and my math is limited. I have different questions: 1. Is it correct to say, that $C: [0,1]^d \rightarrow [0,1]$ is a mapping from a multidimensional distribution to its one dimensional margins? Or that this is a mapping from the unit cube to the uniform margins? 2. I am often reading about the Fréchet-Hoeffding bounds, these seem to basically give the co-domain of the copula? Or what do they state and why do I need them? When are they reached? If I understood it correctly, than the lower bound is reached, when d=2 and the random variables are perfectly negatively dependent and the upper bound is reached, when the random variables are perfectly positively dependent? 3. In every book about copulas I read about the Sklar's Theorem. I do not understand the intuition of this mathematical theorem. What does it say and why do I need it in the context of copulas? 4. The definition of a copula is: It is a multivariate distribution with uniform marginal distributions, but I thought, that I can take any distribution for the marginal distributions and not only the uniform distribution? I know these are beginners questions, but I could not find a good introductory book about copulas. Every book in this topic gives a lot of mathematical definitions and theorems and I am not on this level, unfortunately.


title: Pearson's chi-square statistic for goodness-of-fit in generalized linear models

text: To test the goodness-of-fit in generalized linear models with parameter $\beta$, from a note > Another measure of discrepancy is the generalized Pearson’s chi-square > statistic $$\sum_{i=1}^m \frac{(y_i − \hat{μ_i})^2}{V(\hat{μ_i})}$$ where $\hat{μ_i} := μ_i(\hat{\beta})$, $\mu_i$ is the mean of $Y_i$ under $x_i$, and $V$ is the variance function of the mean. But according to Wikipedia, the Pearson’s chi-square statistic seems to be $$\sum_{i=1}^m \frac{(y_i − \hat{μ_i})^2}{\hat{μ_i}}$$ So why are these two different? Thanks!


title: What do the first k factors from factor analysis maximize?

text: In principal components analysis, the first $k$ principal components are the $k$ _orthogonal_ directions with the maximum variance. In other words, the first principal component is chosen to be the direction of maximum variance, the second principal component is chosen to be the direction orthogonal to the first with the maximum variance, and so on. Is there a similar interpretation for Factor Analysis? For example, I'm thinking that the first $k$ factors are the factors that best explain the _off-diagonal_ components of the original _correlation_ matrix (in the sense of, say, squared error between the original correlation matrix and the correlation matrix defined by the factors). Is this true (or is there something similar that we can say)?


title: Why use ANOVA at all instead of jumping straight into post-hoc or planned comparisons tests?

text: Looking at a between-groups ANOVA situation, what do you get by actually doing such an ANOVA test first, and second do post-hoc (Bonferroni, Šidák, etc.) or planned comparisons tests? Why not skip the ANOVA step entirely? I gather that in such a situation, the one benefit of between-groups ANOVA is to being able to use Tukey's HSD as a post-hoc test. The latter needs the within-groups mean square from the ANOVA table to calculate its related standard error. But, the Bonferroni and Šidák adjustments to unpaired t-tests do not need any ANOVA inputs. I would like to raise the same question regarding within-groups ANOVA situation. I know in such case, the Tukey's HSD test is not a relevant consideration making this question even more pressing.


title: Hypothesis test for summary data involving percentages

text: 500 patients had an operation to remove stomach cancer. There exists a clinical guideline that suggests that a minimum of ten lymph nodes should be removed during the operation. I am therefore interested in the proportion of patients who had at least ten nodes removed. The variable `Fewer_or_more_than_ten` contains this information. I am using R.               # Example data     set.seed(10)     study_data <- data.frame(                   Gender = sample(c("Male","Female"),500,T),                   Fewer_or_more_than_ten = sample(c("Fewer","More"),500,T))      I can tabulate these data to see how many people got an adequate lymph node dissection (i.e., had more than ten nodes removed).               # Tabulate data     table(study_data$Fewer_or_more_than_ten)      I want to know whether this varies with gender, so I can type:               # Stratify by gender     table(study_data$Gender, study_data$Fewer_or_more_than_ten)      What I want to do is to report the percentage of women and men who got an adequate lymph node dissection. So for that I would report the percentage shown under the "More" column in the below table.               # Calculate frequencies     prop.table(table(study_data$Gender, study_data$Fewer_or_more_than_ten), 1) * 100      ### My question is: I would also like to report a p value with these numbers to show whether or not there appears to be a significant difference between the proportion of men versus women who got an adequate dissection. How should I do this? I wasn't sure whether I should try to calculate confidence intervals for the percentages and then compare them to see if they overlap, or if I should apply a Chi-square test or Fisher's exact test to the tabulated data, etc. I am also curious, if gender could take on three values rather than two (e.g., Male, Female, and Unknown), and I wanted to extend my hypothesis test to these three groups, how would that change my approach?


title: How to determine the statistical significance of a score from a general population?

text: Let's say I got a reading of 5 for a given measurement. For a random population, I got the following dataset of readings:               random <- c(2,2,1,1,5,76,3,89,35,66,77,22,99,100,0,34,67,88,95)      I want to determine the probability of getting a reading equal or less than my reading of 5 based on this random population. Could you mind to teach me how to do so? Many thanks.


title: Significance and credibility intervals for interaction term in logistic regression

text: I fitted a Bayesian logistic regression in WinBugs and it has an interaction term. Something like this: $$\mathrm{Prob}(y_{i}=1) = \mathrm{logit}^{-1} (a + b_{1}*x_{i} + b_{2}*w_{i} + b_{3}*x_{i}*w_{i})$$ where $x$ is a standardized continuous variable, and $w$ is a dummy variable. In reality the model is more complicated, but I wanna keep things simple. It happens that the interaction term is "significant", but not the single predictors. For instance, $\mathrm{mean}(b_{1}) = -.2$ and $95%$ quantile: $(-1.3$ and $.7)$ $\mathrm{mean}(b_{2}) = -.4$ and $95%$ quantile: $(-1.3$ and $.5)$ $\mathrm{mean}(b_{3}) = 1.4$ and $95%$ quantile: $(.4$ and $2.5)$ Do you guys have any advice on how to react to this finding? I thought that I could compute 95% credibility intervals for the whole effect of $x$ when $w=1$. This would be: 95% quantile for total effect of x, conditional on $w=1$: $(-1.3+.4$ and $7+2.5) = (-.9 + 3.2)$ Is this correct? If not, what should I do? Any references on the subject? Thanks,


title: Derivation of E step in EM algorithm

text: While im going through the derivation of E step in EM algorithm for pLSA, i came across the following derivation at this page. Could anyone explain me how the following step is derived. $\sum_z q(z) log \frac{P(X|z,\theta)P(z|\theta)}{q(z)} = \sum_z q(z) log \frac{P(z|X,\theta)P(X,\theta)}{q(z)} $


title: How to analyse three independent variables and two dependent variables?

text: I have measured the following variables on 200 individuals.   * Independent variables:     * Measures of personality        * Extraversion       * Conscientiousness     * Leader Member Exchange (LMX) quality   * Dependent variable     * Measures of sales performance        * Sales achieved       * Repeat sales Should I use Multiple regression? I seem to have 3 IVs and 2 DVs, which is causing confusion.


title: What is a "tempered EM algorithm"?

text: In the paper of Probabilistic Latent Semantic Analysis by Hofmann, the author fits the model for _document_ $\times$ _word_ matrix through EM Algorithm in section 3. I was able to follow the derivation and meaning of the model derived in it. However in the later section, the author mentioned about Tempered EM for improving the generalization. Could anyone explain or point me to the location where I can understand the actual meaning of Tempered EM Algorithm.


title: Predicting forecasts for next 12 months using Box-Jenkins

text: I am building a Box-Jenkins model in Excel using solver. The model is AR(2). The data that I have contains trend and seasonality both. I know how to remove seasonality using seasonal indexes and add it back to the forecast. But, how do I handle trend? If I remove trend from the data, how should I add it back to the forecast? Also, is the excel solver best way to find the AR parameters?


title: Tobit model with t-distribution

text: I've got a question concerning the estimation of a tobit model with the AER package in R. I observed t-distributed residuals, so that the assumption of normal distributed std errors is violated. Fortunately it's possible to choose "t" as the distribution when fitting with `tobit()`. It is `tobit(<formula>, dist="t")`. After changing the distribution from t to normal, the AIC increases, so this seems to be a good idea. My problem is: Is it possible to flip the normal distribution that way? Is it still a tobit model and do I have to care about something when changing the distribution in a tobit model? I appreciate any thought!


title: Constructing a tridiagonal Markov Chain given a desired distribution

text: Given a desired distribution over a set of finite states, I wish to construct a tridiagonal Markov chain that would have a stationary distribution same as the desired distribution. The reason I want to do this is as follows. I have a linear sequence of blocks that an agent can traverse. Each time, the agent can choose to stay at the current block, move forward one block, or move backward one block. However, I also have a target distribution such that over time, I want the frequency that each block gets visited to be similar to the target distribution. I personally don't know any existing methods/tools that could easily solve this due to my limited knowledge in statistics. I'd appreciate it if anyone can provide some advice on how I might be able to approach this problem.


title: Training a hidden Markov model

text: When it comes to the training of a Hidden Markov Model using multiple training instances should I first take a single instance and train the model until the convergence and then move on to the next instance or should I use multiple instances sequentially to train the model?


title: Does adjusted R-square seek to estimate fixed score or random score population r-squared?

text: Population r-square $\rho^2$ can be defined assuming fixed scores or random scores:   * Fixed scores: The sample size and the particular values of the predictors are held fixed. Thus, $\rho^2_f$ is the proportion of variance explained in the outcome by the population regression equation when the predictor values are held constant.   * Random scores: The particular values of the predictors are drawn from a distribution. Thus, $\rho^2_r$ refers to the proportion of variance explained in the outcome in the population where the predictor values correspond to the population distribution of the predictors.  I've previously asked about whether this distinction makes much difference to estimates of $\rho^2$. I've also asked generally about how to calculate an unbiased estimate of $\rho^2$. I can see that as the sample size gets larger the distinction between fixed- score and random-score gets less important. However, I'm trying to confirm whether adjusted $R^2$ is designed to estimate fixed score or random score $\rho^2$. ### Questions   * **Is adjusted $R^2$ designed to estimate fixed score or random score $\rho^2$?**   * **Is there a principled explanation of how the formula for adjusted r-square relates to one or other form of $\rho^2$?** ### Background to my confusion When I read Yin and Fan (2001, p.206) they write: > One of the basic assumptions of the multiple regression model is that the > values of the independent variables are known constants and are fixed by the > researcher before the experiment. Only the dependent variable is free to > vary from sample to sample. That regression model is called the _fixed > linear regression model_. > > However, in social and behavioral sciences, the values of independent > variables are rarely fixed by the researchers and are also subject to random > errors. Therefore, a second regression model for applications has been > suggested, in which both dependent and independent variables are allowed to > vary (Binder, 1959; Park & Dudycha, 1974). That model is called the random > model (or correction model). Although the maximum likelihood estimates of > the regression coefficients obtained from the random and fixed models are > the same under normality assumptions, their distributions are very > different. The random model is so complex that more research is needed > before it can be accepted in place of the commonly used fixed linear > regression model. Therefore, the fixed model is usually applied, even when > the assumptions are not met completely (Claudy, 1978). Such applications of > the fixed regression model with assumptions violated would cause > "overfitting," because the random error introduced from the less-than- > perfect sample data tends to be capitalized in the process. As a result, the > sample multiple correlation coefficient obtained that way tends to > overestimate the true population multiple correlation (Claudy, 1978; Cohen & > Cohen, 1983; Cummings, 1982). So I wasn't clear whether the above statement is saying that adjusted $R^2$ compensates for error introduced by the random model or whether this was just a caveat in the paper flagging the existence of the random model, but that the paper was going to focus on the fixed model. ### References   * Yin, P., & Fan, X. (2001). Estimating $R^2$ shrinkage in multiple regression: A comparison of different analytical methods. The Journal of Experimental Education, 69(2), 203-224. PDF


title: Is any quantitative property of the population a "parameter"?

text: I'm relatively familiar with the distinction between the terms statistic and parameter. I see a statistic as the value obtained from applying a function to the sample data. However, most examples of parameters relate to defining a parametric distribution. A common example is the mean and standard deviation to parameterise the normal distribution or the coefficients and error variance to parameterise a linear regression. However, there are many other values of the population distribution that are less prototypical (e.g., minimum, maximum, r-square in multiple regression, the .25 quantile, median, the number of predictors with non-zero coefficients, skewness, the number of correlations in a correlation matrix greater than .3, etc.). Thus, **my questions** are:   * Should any quantitative property of a population be labelled a "parameter"?   * If yes, then why?   * If no, what characteristics should not be labelled a parameter? What should they be labelled? And why? ### Elaboration on confusion The Wikipedia article on estimators states: > An "estimator" or "point estimate" is a statistic (that is, a function of > the data) that is used to infer the value of an unknown parameter in a > statistical model. But I can define the unknown value as .25 quantile and I can develop an estimator for that unknown. I.e., not all quantitative properties of a population are parameters in the same way that say the mean and sd are parameters of a normal distribution, yet it is legitimate to seek to estimate any quantitative population property.


title: How to compute multiple date ranges for wildlife occurrences

text: Say we are trying to predict occurrences of lightning bugs based only on the dates on which people reported seeing them. So climate change, weather, and other influencing factors are imaginary. So there are $n$ recorded occurrences of lightning bugs. Since many people have contributed and the recording is not standardized, likely only the days of the occurrences - and not the number - should be considered. If every day on which a lightning bug was spotted was within July and August, and are densest around July 20, then it seems simplest to return the first and second quartile. But say there was a high density of occurrences in May, and then none again until July and August. The quartile method would be false. Where should I start in writing an algorithm to return multiple date ranges?.   Could a probability be assigned to the result?.


title: statistical test to see if relationship is linear or non-linear

text: I have an example data set as follows:               Volume <- seq(1,20,0.1)     var1 <- 100      x2 <- 1000000     x3 <- 30          x4 = sqrt(x2/pi)     H = x3 - Volume     r = (x4*H)/(H + Volume)          Power = (var1*x2)/(100*(pi*Volume/3)*(x4*x4 + x4*r + r*r))          Power <- jitter(Power, factor = 1, amount = 0.1)     plot(Volume,Power)      ![enter image description here](http://i.stack.imgur.com/vHBcY.png) From the figure, it can be suggested that between a certain range of 'Volume' and 'Power' the relationship is linear, then when 'Volume' becomes relatively small the relationship becomes non-linear. Is there a statistical test for illustrating this? With regards to some of the recommendations shown in the responses to the OP: The example shown here is simply an example, the dataset I have looks similar to the relationship seen here although noisier. The analysis that I have conducted so far shows that when I analyse a volume of a specific liquid, the power of a signal drastically increases when there is a low volume. So, say is I only had an environment where the volume was between 15 and 20, it would almost look like a linear relationship. However, by increasing the range of points i.e. having smaller volumes, we see that the relationship is not linear at all. I am now looking for some statistical advice on how to statistically show this. Hope this makes sense.


title: Study design retrospectively comparing detection rates of two diagnostic tests using case-control-like design

text: I am trying to find a name for the kind of study design/analysis I am working on so I can use it to look up best practices and examples of people who have done this kind of thing well. Essentially I am trying to retrospectively compare two ways of taking biopsies of an organ. Since we are taking biopsies, we assume that the PPV of each test is 100%, so the question is simply "does biopsy method A find more cancer than biopsy method B?" We have lots of examples of biopsy method B but only a few of A, so we would like to take patients who got biopsy type A and match them to one or more patients who got biopsy type B based on basic demographic information and risk factors for the disease before doing the comparison. This sounds a bit like a case-control study, but it's not quite. Is there a name for what I am trying to do here? Can anyone point me to a paper that talks about or uses this kind of design?


title: Model selection in nlme's

text: I can think of four ways to perform model selection `nlme`'s: `LRT`, `AIC`, `BIC` and `CV`. It would be much appreciated if anyone could add to/correct this list. I vaguely recall somebody computing `cAIC` in `SAS`, but not how df was defined?


title: Resources for matrix calculus for optimization

text: I'm a grad student trying to absorb from the book Pattern Recognition and Machine Learning. However, I found that I really need a good grasp of matrix calculus before I can deduct the formulas myself (since I think in this way learning could be more effective). For example, when deducting Gaussian Mixture Model using EM algorithm, I can't really do the M-step on my own since I don't know how to do calculus on vectors/matrix. What I wish to solve is to be able to solve derivatives like $x^TAx$ with respect to $x$. I did read wikipedia and know the basics of the idea of taking derivatives with respect to vector, but I hope to get a sense and make myself more confirmable with these operations. Questions: 1) Does people memorize the derivatives of the basic formulas (like $x^TAx$, or derivaitive of covariance) and then expand?   2) What good resources especially books with practices problems could be recommended? Thanks a lot


title: How does one derive the Fisher-consistency correction factor for Rousseeuw and Croux's $S_n$ (and $Q_n$)?

text: In "Alternatives to the Median Absolute Deviation" (Rousseeuw and Croux, _J. Amer. Statistical Assoc_ , **88** (424), 1993, pp.1273-1283) and a few other papers from the same authors published in 1992-1993, Rousseeuw and Croux introduce their robust scale metrics $S_n$ and $Q_n$. They offer correction factors for finite sample size and for Fisher consistency for the normal distribution. In particular, they find that the correction factor $c$ for Fisher consistency for the normal distribution satisfies $$\Phi\left(\Phi^{-1}(3/4)+c^{-1}\right)-\Phi\left(\Phi^{-1}(3/4)-c^{-1}\right) = \frac{1}{2}.\tag{1}$$ In words, $c$ is such that the interval centered at $\Phi^{-1}(3/4)$ and with width $2 c^{-1}$ covers exactly 50% of the probability mass of the normal distribution. On the other hand, more or less by definition of $S_n$, we have (for any distribution) $$c^{-1} = \text{med}_X (\text{med}_Y \ \left\vert X - Y \right\vert),\tag{2}$$ where $\text{med}_Z$ stands for the median over some random variable $Z$. In simulation experiments, I found $(1)$ to be true, but I don't see how the argument works, and the paper doesn't offer any help. Does anyone see the trick? My hope is that the argument generalizes to other distributions sometimes. Indeed it seems (from experiments) that for $X \sim \text{Exp}(\lambda)$, we have $$F\left(F^{-1}(r)+c^{-1}\right)-F\left(F^{-1}(r)-c^{-1}\right) = \frac{1}{2}$$ with $r \approx 0.6$ (but for the beta distribution there is no such formula). Is it known if this is true exactly for any $r$, and if so, what is that $r$? More generally, are there other ways than equation $(2)$ to obtain a formula (possibly implicit, such as $(1)$) for $c$ for non-normal distributions?


title: How to compute precision/recall for multiclass-multilabel classification?

text: I'm wondering how to calculate precision and recall measures for multiclass multilabel classification, i.e. classification where there are more than two labels, and where each instance can have multiple labels?


title: Artificial dataset generator for classification data

text: I would like to generate some artificial data to evaluate an algorithm for classification (the algorithm induces a model that predicts posterior probabilities). These are some basic properties of the dataset:   * Features have to be continuous   * Response variable is dichotomous (either 0 or 1) I would like to test whether the algorithm can cope with:   * Many feature / high dimensional problems   * noise (it can drop features)   * Multi-modality    * ??? (how do I simulate correlation etc.) I intend to implement the algorithm in R or Matlab. I can sample from multivariate normal distributions and specify a covariance matrix. I would appreciate any feedback.


title: Which R package can estimate parameters in my Kalman Filter

text: I am using Kalman Filter to estimate my state variables for a time series data. But all my parameters in KF are unknown.               S(t) = A*S(t-1) + B*X(t) + Q   ----------Sate process         Z(t) = H*S(t)   + W*D(t) + R   ----------Measurement process where,           T: the total number of time slots in the time series data         N: number of state variable         M: number of measurement variable         S: state variable, N*T dimension matrix         Z: measurement variable, M*T dimension matrix         Q: state process noise         R: measurement process noise         X: input in state process      L dimension for a time slot         D: input in measurement process, I dimension for a time slot         A, B and H: parameter matrix       I want to solve                1. Best A, B, and H matrix for the KF. All parameters in the model are constrained.      2. Estimate all state variables after getting best A, B and H.      One example I want to solve is to estimate sales effect from different media channels and media synergy effects.                Z(t) is the sales number for each week       S(t) is the sales effect from each media channle and media synergy effect; a sales effect from a media channel will be a state variable.      X(t) is the media spending      D(t) are seasonality, macro-economy factor and other extra independent varialbes to the model      I have checked some R packages. It seems that some of them have specific KF formulation. Which R packages can I use to solve the two problems listed above? Do I need to call optim() together with those R packages in order to solve problem 1? Thanks!


title: DLM package, state equation, maximum likelihood with constant terms

text: My question is about the DLM package and the dlmMLE. Lest's say that I have a bivariate model of this kind: $$Y(t)= Fb(t)+e(t)$$ $$b(t) = u + Gb(t) +w(t)$$ $$b(t)=(b_1(t),b_2(t)),\quad u = (u_1,u_2)$$ Rewriting the model in matrix notation, the state equation becomes: $$d(t) = H d(t) +w(t),$$ where $d(t)= (b(t),1,1)$ and $H$ is a matrix 4x4 given by the merging the 2 matrices 2x4 $[G,\text{diag}(u)]$ and $ [0,\text{diag}(1,p)]$, and the variance of the innovations $w$ is a matrix of zeros except for the upper right 2x2. Morover I assume that $C_0=\text{diag}(1e3,1e3,0,0)$ Since it doesn't work well I was wondering if it is the right way to write it? Thanks I would appreciate any suggestions.


title: An investment and variance question for monthly payments

text: I have a question regarding a financial/statistical problem. How do you calculate the variance of the outcome of an investment in a stock, when the investment is so called time diversified, i.e. spread out over time. So, an investor invests an amount every month during 24 months in the same stock. What is the variance of the outcome after 24 months (normal distribution assumed)? Expected yearly return 5%, annual std 10%. I have spent a few hours googling but have not found anything (at least I could) translate into solving this. Thanks in advance Toby


title: Multinomial logistic regression - Behaviour biology

text: I am a student in biology, currently finishing my master in Behavior, evolution and conservation. I have been in the Swiss national park for my field work and I have some trouble to analyze my data. So basically what I have is a lot of factors and what I need to do is a multinomial logistic regression : Factors :   * Behaviour (4 states - Moving, Feeding, Resting, Runing)   * Age (Of the individual)   * Temperature   * Valley (where the chamois was, 2 different choices)   * Year (of the observation)   * Month (of the observation)   * Kid (if it has a kid or not)   * Individual (The number written on the tag he had on his ear) What I want to do is to check what factors influence the Behavior. I tried several things with R but can't use the `mlogit` package which is the one I have been told to use. I have also tried in JMP, now the problem I have here is that I just clicked on "Fit Y in function of X" and selected my response variable and my factors and bam, I had results, but to be honest this seems very simple and I was wondering if I am not missing something somewhere. Edit : Here is what R returns when I try the `mlogit` function :               mlogit(Behavior~Age+Temp+Valley+Individual+Year+Month, Merge)     Error in `row.names<-.data.frame`(`*tmp*`, value = value) :        invalid 'row.names' length      Can anyone here that can help me either with R or with JMP?


title: cor() standard deviation error is zero using random sample of data

text: I'm reading a paper by Max Kuhn called Building Predictive Models in R Using the caret Package. He uses two data sets, `desc` and `mutagen`. Their content is not important here. He runs the following code:               library("caret")     set.seed(1)     inTrain <- createDataPartition(mutagen, p=3/4, list=FALSE)     trainDescr <- descr[inTrain,]     descrCorr <- cor(trainDescr)     highCorr <- findCorrelation(descrCorr, 0.90)     trainDescr <- tranDescr[, -highCorr]      When I try to run `cor(trainDescr)`, I get the error               In cor(trainDescr) : the standard deviation is zero      followed with               Error in if (mean(x[i, -i]) > mean(x[-j, j])) { :        missing value where TRUE/FALSE needed      when running `findCorrelation()`. If I try and run the code on the whole data set `descr`, and not the random sample, everything works fine. I'm running R version 4.12.1 on 32-bit ubuntu. Hope this makes sense to someone.


title: Logistic Regression Cost Function issue in Matlab

text: I'm trying to implement a logistic regression function in matlab. I calculated the theta values, linear regression cost function is converging and then I use those parameters in logistic regression function as a decision boundary. I'm trying to obtain an overfit logistic regression tree to show how cost function behaves during overfitting with respect to training set size. I get a converging training error which is expected. I also get an increasing validation error which is also expected but even though they are increasing, my validation error values are always below the training error values. Is it suppossed to be like that (because we are using the logarithm of the hypothesis function and not the squared error function) or am I doing something wrong?. (I can't post images because I'm a new member) My cost function calculations are included in the code part, I'd appreciate any help. Thanks!                   for i=1:row         if y(i)>0             classification(i)=1;         else             classification(i)=0;         end     classpredict(i)=1/(1+exp(-(ypredict(i))));     end          for i=1:size(data,1)-row         if validy(i)>0             classvalid(i)=1;         else             classvalid(i)=0;         end     classvalidpredict(i)=1/(1+exp(-(validpredicty(i))));     end               for i=1:size(classpredict)                  costclasstrain=-(costclasstrain+(classification(i))*(log(classpredict(i)))+(1-classification(i))*(log(1-classpredict(i))))/row;     end          for i=1:size(classvalidpredict)                costclassvalid=-(costclassvalid+(classvalid(i))*(log(classvalidpredict(i)))+(1-classvalid(i))*(log(1-classvalidpredict(i))))/(size(data,1)-row);     end


title: Difference between survdiff log-rank and coxph log-rank

text: I'm using the survival package in R to analyze clinical data. I am analyzing two different groups of patients, when I calculate survdiff in order to compare the curves, I got p= 0.135, but when I adjust the model using coxph and different covariates, let say clinical cancer stages, , I got an overall logrank score of 0.0005793 for 5 covariates. My question is, could I use this late logrank p-value to say that adjusting the model with more covariates the difference between the curves is signifficative? here is the data               survdiff(formula = my.surv ~ final_table$G)          n=56, 14 observations deleted due to missingness.                           N Observed Expected (O-E)^2/E (O-E)^2/V     final_table$G=1  4        2     1.43    0.2294     0.247     final_table$G=2 52       24    24.57    0.0133     0.247           Chisq= 0.2  on 1 degrees of freedom, **p= 0.619**       And this is the coxph results               coxph(formula = Surv(final_table$Time_surv, final_table$Survival) ~ final_table$G + final_table$ST)            n= 56, number of events= 26        (14 observations deleted due to missingness)                             coef exp(coef)  se(coef)     z Pr(>|z|)     final_table$G2    2.094e-01 1.233e+00 7.532e-01 0.278    0.781     final_table$STII  1.883e+01 1.501e+08 5.739e+03 0.003    0.997     final_table$STIII 1.998e+01 4.773e+08 5.739e+03 0.003    0.997     final_table$STIV  2.089e+01 1.186e+09 5.739e+03 0.004    0.997                        exp(coef) exp(-coef) lower .95 upper .95     final_table$G2    1.233e+00  8.111e-01    0.2817     5.396     final_table$STII  1.501e+08  6.662e-09    0.0000       Inf     final_table$STIII 4.773e+08  2.095e-09    0.0000       Inf     final_table$STIV  1.186e+09  8.430e-10    0.0000       Inf          Concordance= 0.74  (se = 0.057 )     Rsquare= 0.37   (max possible= 0.957 )     Likelihood ratio test= 25.86  on 4 df,   p=3.381e-05     Wald test            = 4.02  on 4 df,   p=0.4033     Score (logrank) test = 19.67  on 4 df,   **p=0.0005793**      Thanks Thanks to comments I did this analysis, survdiff with roup and stage               survdiff(formula = Surv(final_table$Time_surv, final_table$Survival) ~      final_table$G + final_table$ST)          n=56, 14 observations deleted due to missingness.                                           N Observed Expected (O-E)^2/E (O-E)^2/V     final_table$G=1, final_table$ST=III  3        2    1.149     0.630     0.668     final_table$G=1, final_table$ST=IV   1        0    0.279     0.279     0.285     final_table$G=2, final_table$ST=I   15        0    8.715     8.715    13.547     final_table$G=2, final_table$ST=II   2        1    1.816     0.367     0.402     final_table$G=2, final_table$ST=III 30       19   13.067     2.693     5.540     final_table$G=2, final_table$ST=IV   5        4    0.973     9.413     9.935          Chisq= 23.2  on 5 degrees of freedom, p= 0.000313      So the final value is totally significant, but now I got 6 curves, more or less this is what I want, how the group and the stage is affecting the survival. What do you think?


title: Estimating latent performance potential based on a sequence of observations

text: **Context**   * you have 200 observations of an individual's running time for the 100 metres measured once a day for 200 days.   * Assume the individual was not a runner before commencement of practice   * Based on the observed data and the 199 other observations, you want to estimate the latent time it would take the individual to run if they (a) applied maximal effort; and (b) had a reasonably good run for them (i.e., no major problems with the run; but still a typical run). Let's call this **latent potential**. Of course, the actual data would not measure latent potential directly. The data would be noisy:   * Times would vary from run to run   * On some days the individual would be particularly slow because of one or more possible problems (e.g., tripping at the start, getting a cramp half way through, not putting in much effort). Such problems would result in massive outliers   * On some days the individual would be slower than you'd expect, perhaps because of more minor issues.   * In general, with practice the runner would be expected to get faster in latent potential.   * In rare cases, it is possible for the runner to get slower in latent potential (e.g., injury) The implications of this:   * The occasional slow time might provide minimal information on what the individual is capable of.    * A fast time for the individual suggests that the individual is capable of such a fast time, but a small amount of this fast time might be good fortune on the day (e.g., the right wind, a little luck on the start). **The question:** Thus, how could one estimate latent potential at each of the 200 time points based on the available data and a few assumptions about the nature of running times? **Initial Thoughts:** I imagine there would be some form of Bayesian approach that combined the available information and assumptions to form an estimate, but I'm not sure where to look for such models. I'm also not quite clear how the effectiveness of such a model would be evaluated.