title: Interpretation of plotted component smooth functions of generalized additive model - transformed variables

text: I have a plot of a GAM which models deforestation severity. ![Component smooth functions of a GAM](http://i.stack.imgur.com/A9tXG.png) The explanatory variables are transformed, with log10, inverse, and square root. Is it meaningful to interpret a gradient like this: 'Deforestation severity increases as log10(Distance to forest edge) increases' or 'Deforestation severity increases as distance to forest edge increases' for plot e (note, the lower the value of the y axis the greater the deforestation severity). If not, how could I plot the component smooth functions untransformed?


title: Recursive Feature Elimination Fails to Output as Expected

text: Currently I am using rfe function in the "caret" package to do feature selection. There are 380 variables as input candidates. I have done many trials and I noticed that something weird always happens, as the rfe tends to output all candidate variables as the best subset while I have set the parameter size to be 5 to 10. I am wondering whether anyone else has come across this problem. Any help will be appreciated. Thanks in advance for any help.


title: Ratio of lengths of two confidence intervalls

text: I have two random variables: (1) With standard normal distribution. Confidence interval $I_1$, which is centered and has probability of $(1-\alpha)$ (2) With T-distribution. Confidence intervall intervall $I_2$, also centered and has prob. $(1-\alpha)$. First I had to compute the length of both: I got $\mathcal{L}(I_1) = 2 \cdot z_{1-\alpha/2} \frac{\sigma}{\sqrt{n}}$ and $\mathcal{L}(I_2) = 2 \cdot t_{1-\alpha/2;n-1} \frac{S}{\sqrt{n}}$. Hope this is correct. Now I have to calculate the asymptotic distribution of the ratio: $\Lambda = \frac{\mathcal{L}(I_2)^2}{\mathcal{L}(I_1)^2}$. As a hint I got $(n-1)\frac{S^2}{\sigma^2} \sim \chi_{n-1}^2$. So far I got: $\Lambda = \frac{z_{1-\alpha/2}^2}{t_{1-\alpha/2;n-1}^2}\cdot\frac{S^2}{\sigma^2}$. How do I express the ratio of the two distributions? I know that t-Distribution gets to a standard normal if n goes to infinity, but does that help? And do I get the $(n-1)$ from the hint somehow, or is this not really relevant?


title: A resource on concepts underlying statistics, not the techniques used in applied stats

text: What book is the most thorough treatment of fundamental concepts in statistics? I am not asking for a book on details of the methods of calculations and procedures, I am mainly interested in a book that thoroughly explains the foundational concepts ... an intuitive/illustrated/visual approach to the core ideas ... rather than load of mathematics equations etc. Size of the book is no problem ... even a large multi-volume text would do ... heck even a web resource would do.


title: Handling stationarity issues in proc ucm/state space time series models

text: Hope I'm able to find someone who can answer this question. The previous one didn't get answered! Proc ucm is the SAS implementation (using state space concepts) to isolate the unobserved trend, seasonality & estimate the coefficients of regressors simultaneously. Documentation/research on proc ucm is sparse. There are two questions I'm trying to find answers for: a. Does the dependent variable have to be treated for weak stationarity (differenced) before running proc UCM? Ans. http://www.iasri.res.in/sscnars/socialsci/17-Structural%20Time%20Series%20Models%20for%20Describing%20Trend%20in%20All%20India%20Sunflower%20Yield%20Using%20SAS.pdf The attached paper seems to suggest that I don't need to worry about stationarity while using proc ucm. However, I have found a contrary situation in my data. Consider my dependent is y. And, t & s are the smoothed trend & season components isolated by proc ucm. I'd have expected the series (y-t-s) to be stationary. But, it was not. Thereby, I conclude that proc ucm is not capable of handling non-stationary time series, until I explicitly difference/stationarize the dependent. Is this right? b. I also have regressors which exhibit co-integrated relationships with y. From Granger's research paper, it is evident that spurious regression results when co-integrated series exists. But, I hear that's only if we use OLS-based proc reg. Proc UCM is based on maximum likelihood estimation. Is co- integration not a problem with maximum likelihood estimation based regression? I could be ambiguous in my problem statement above and can clarify if need be.


title: Model R&D decision and intensity

text: What is the preferred model to analyze firms' R&D decision and intensity?   Normally the questionaires include two steps:   1. firms are asked whether they engage in R&D (binary variable (yes/no))   (around 70 % answer no)   2. If yes, firms report how much they invest. Hence, 70% of all firms report 0 investment, while 30% > 0. I have seen several papers using different models.   1. Simple OLS    2. Heckman selection model (normally without an exclusion restriction)   3. Zero-inflated count model (Possion or negative binominal)   4. Hurdle model   5. ? What is the best model for this kind of data?   How do I know that this particular model is the most appropriate?


title: Sequential Inference And Evidence (Jaynes 2003): Is it valid? Is it used?

text: Exploring the work of ET Jaynes, _Probability Theory_ (11th Printing 2013) has led to consideration of the technique he identifies as **Sequential Inference** (p. 96); where the _evidence_ , in decibels, accumulates until the investigator either (1) stops with acceptance, (2) stops with rejection or (3) continues with another test. Jaynes' technique appears to be a special case of, and perhaps a more Bayesian approach than, sequential analysis as proposed by Abraham Wald. But I am certainly not an expert here. This looks to be an incredibly powerful technique that can be used to investigate a claim of compliance/non-compliance in a cost effective and time efficient manner. Yet when searching the internet and also checking other Bayesian texts (such as Silvia and Skilling) there is little depth on it, if it is mentioned at all. So these questions are posed: (1) Is it used in practice? (2) If so, in what way, if not, why not? (3) Are there critical issues/pitfalls with it's use in practice? (4) Is (are) there any in-depth reference(s) with case studies of application? **Addendum** Since posting this earlier, we have found an excellent collection of information here. Coming from an engineering background, by nature we tend to be data driven, and this would put Engineers solidly in the Frequentist camp. However, time and again we have seen that these analyses should start from a position of logic, and this is what a _proper_ Bayesian approach will accomplish, it is clear. So there is no doubt in our mind that the Bayesian method _trumps_ the Frequentist one. However, what we are specifically interested in is the concept of _evidence_ in _decibels_ (it seems that _bels_ or maybe the _Jaynes Scale_ (0,10] analogous to the _Richter Scale_ might have been a more interesting way to consider it). But is this approach used for more routine issues found in practice, and are there any case studies of such?


title: conditional logit model will not converge

text: I am running a conditional logit model and it will not converge! The DV has 6 levels of vote choice, essentially parties one can vote for, and I've got one predictor, individual level economic perception. Because the independent variable does not vary across alternatives, I created interaction terms by multiplying the options to the economic perceptions. There has to be something wrong in the code, because I keep getting the error message that it "ran out of iterations and will not converge". I thought I could take out one of the constituent terms, but the does not fix things either... Any help would be greatly appreciated :) The code I used is below:               choice<- as.numeric(data1$voteChoice1)     table(choice)     selection <- c()     econ <- c()     chooser <- c()          c <- 1          for(i in 1:length(data1$voteChoice1)){           selection[c:(c+5)] <- as.numeric(data1$voteChoice1)[i]           econ[c:(c+5)] <- data1$economy[i]       chooser[c:(c+5)] <- i       c <- c + 6     }          econ <- data1$economy     vote <- choice == selection          otherecon <- na.omit(econ* as.numeric(choice=='other'))     grnecon<- na.omit(econ * as.numeric(choice=='green'))     cduecon <- na.omit(econ * as.numeric(choice=='cdu'))     fdpecon<- na.omit(econ * as.numeric(choice=='fdp'))     spdecon<- na.omit(econ * as.numeric(choice=='spd'))     dlecon<- na.omit(econ * as.numeric(choice=='dl'))          other <- na.omit(as.numeric(choice=='other'))     grn <- na.omit(as.numeric(choice=='green'))     cdu <- na.omit(as.numeric(choice=='cdu'))     fdp<-na.omit(as.numeric(choice=='fdp'))     spd <- na.omit(as.numeric(choice=='spd'))     dl <- na.omit(as.numeric(choice=='dl'))     votechoice <- na.omit(data1$voteChoice1)          modelc <- clogit(vote ~ grn + cdu + other + fdp + spd + dl + grnecon + cduecon + fdpecon + spdecon + dlecon  + otherecon + strata(chooser), na.action=na.omit)     summary(modelc)


title: How to choose data for training a predictive model for attrition prediction

text: Trying to build a predictive model for attrition prediction at service desk/call center. Have daily data on the following parameters: 1.Call quality - QTM (0-100%), 2.No. of calls - Calls(Number) 3.Attendance 4.Customer feedback(1/0) Q1,Q2 (0-100%) for both, agents who left the job and for the ones who are still there, for a duration of 6 months. Aim: to predict agents tendency/probability of staying/leaving based on his/her daily performance. Doubts i have, 1\. how should i use the data to train the model(logistic regression) should it be trained based on the avg of the parameters taken over a duration of 6 months. **if so can we test the daily metrics based on a model which is trained using mean of the parameters for 6 months. Please advise. this is my first attempt at making a predictive model,i have gone thru various case studies/models such as the titanic survival model using logistic regression,Wisconsin DEWS model. I decided to model using the weekly aggregates of the the two populations(attrites and Non-attrites). The Data Set (approx 5 months data,with weekly aggregates of the two populations i.e Attrites and Non-Attrites.) AW1 : Week1 Aggregates of the performance metrics for Attrites NAW1: Week1 Aggregates of the performance metrics for Non-Attrites ![Data-set](http://i.stack.imgur.com/EUMCI.png) Post this i ran a logistic Regression on 80% of this data-set and kept aside the other 20% for testing. Results of the logistic regression: ![results of the logistic regression model.](http://i.stack.imgur.com/wyB8B.png) and then i used the predict function on the 20% of the data which contained 3 data points for both attrites and Non-attrites,so to be 100% accurate the model should have predicted 3 as attrites and 3 as Non-attrites but the correct prediction is 5/6 that is one wrong prediction out of 6. Please help me in interpreting the meaning of the results of the model all the z values are zero im not sure what that signifies. Googled a little regarding the z values = 0 issue and came across some posts on stackoveflow that suggested using "bayesglm" instead of "glm" did that and the results are good at the first look but being a newbie in the field i would like you to guide me with respect to the statistical significance of the issue and is the model really as good as the results of the "bayesglm" or is it just by fluke. ![results with bayesglm](http://i.stack.imgur.com/Y2Cr8.png) the model gives a 100% accurate prediction on the test set now 6/6.


title: Generalized Linear Mixed Model in R with repeated measures

text: I am trying to investigate how four variables (var1=continuous, var2=factor, var3=factor, var4=continuous) influence the number of trials individuals approached (out of total nr of trials --> binomial) across two conditions that differed in food availability (food availability 1 = 42 trials; food availability 8 = 35 trials) (n = 19 individuals). The response variable is binomial as it is the number of trials out of total number of trials. I am using the 'lmer' function of the lme4 package. I thought the additive model I should run would be with random factor ID:               glmer(cbind(appr_Y,appr_N) ~ Condition+Var1+Var2+Var3+Var4+(1|ID), data=dataset,            family=binomial)      However, the result I get shows that Condition is super significant (p < 2e-16) while the other variables aren't, while exploring the data visually shows no difference in the response variable for Condition and the variables having strong effects. Below a dummy representing the large data table:               Con ID  Var1  Var2  appr_Y  appr_N  Trial_total     1   1   10      y   14      6       20     1   2   4       y   10      10      20     1   3   5       n   5       15      20     1   4   32      n   18      2       20     1   5   11      y   3       17      20     2   1   10      y   20      5       25     2   2   4       y   10      15      25     2   3   5       n   24      1       25     2   4   32      n   11      14      25       2   5   11      y   7       18      25      What am I doing wrong? **update** : I analysed the data with GenStat (which doesn't show AIC values) and the output is totally different. In GenStat it asks for the random factor (here ID) and the denominator (here Trial_total), which is different than putting in Appr_Y, Appr_N. **update2** : The above dataset was just a dummy. I hereby provide the 'summary' of the model and the information about the dataset:               > summary(GLMM1)     Generalized linear mixed model fit by the Laplace approximation      Formula: cbind(Appr_Y, Appr_N) ~ Condition + Var1 + Var2 + Var3 + Var4 + (1 | ID)      Data: dataset      AIC   BIC logLik deviance     102.1 113.5 -44.04    88.08     Random effects:     Groups Name        Variance Std.Dev.     ID     (Intercept) 0.59495  0.77133      Number of obs: 38, groups: ID, 19          Fixed effects:                       Estimate Std. Error z value Pr(>|z|)         (Intercept)       -2.43536    0.60237  -4.043 5.28e-05 ***     Condition8         1.14942    0.12274   9.365  < 2e-16 ***     Var1               0.04524    0.04002   1.130   0.2583         Var2Paired        -0.35299    0.47970  -0.736   0.4618         Var3no             0.55914    0.44095   1.268   0.2048         Var4               0.11996    0.06282   1.909   0.0562 .       ---     Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1           Correlation of Fixed Effects:                 (Intr) Cndt8- Var1 Var2P Var3no     Cndtn8-strn -0.128                                 Var1        -0.294  0.015                          Var2unp     -0.474 -0.015 -0.352                   Var3no      -0.178  0.016 -0.310 -0.097            Var4        -0.664  0.021 -0.078  0.467 -0.134     > str(dataset)     'data.frame':   38 obs. of  9 variables:     $ ID          : Factor w/ 19 levels "39","40","41",..: 1 2 3 4 5 6 7 8 9 10 ...        $ Appr_Y      : num  3 12 0 7 27 6 12 1 5 17 ...     $ Appr_N      : num  39 30 42 35 15 36 30 41 37 25 ...         $ Var2        : Factor w/ 2 levels "paired","unpaired": 2 2 2 2 1 1 2 1 2 1 ...     $ Var1        : num  2 16 19 18 13 11 14 1 8 9 ...         $ Var3        : Factor w/ 2 levels "yes","no": 2 2 2 1 2 2 2 1 1 2 ...     $ Var4        : num  2.6 6.87 2.4 1.1 4.32 ...         $ Condition   : Factor w/ 2 levels "1","8": 1 1 1 1 1 1 1 1 1 1 ...     $ n           : num  42 42 42 42 42 42 42 42 42 42 ...      Do I perhaps have to do something with weighing the data as trial nr is not the same across conditions? Or using Appr_Y, total nr of trials instead of Appr_Y, Appr_N ?


title: Variable significance for assessment results

text: I am trying to analyze multiple variables of our student population as it comes to assessment outcomes. I have the assessment change score and demographic information such as (race, sex, age, program, # of programs, time of service, # of services). I am trying to ascertain which variables are significant to the change scores, and how significant.


title: Logistic Regression : How to obtain a saturated model

text: I just read about the deviance measure for the logistic regression. However, the part that is called saturated model is not clear to me. I did an extensive Google search but none of the results answered my question. So far I found out that a saturated model has a parameter for each observation which as a consequence results in a perfect fit. This is clear to me. But: further the fitted values (of a saturated model) are equal to the observed values. Since from my knowledge, logistic regression is used for classification the given observed data are covariates with additional labels $y \in \\{0,1\\}$. However, the deviance measure employs probabilities but not the actual labels. One applies the calculated predicted probability of the logistic regression versus observed probabilities. However, since one has only given labels instead of probabilities I'm confused how to construct a saturated model from these labels?


title: How to aggregate/deal with time series data in multiple markets?

text: I have our company's daily revenue data in multiple markets. The essence of the project is to see how different markets compare to each other based on a few variables, namely the age of the market and its population. I've created 30 day rolling means for each market and plotted them, ![enter image description here](http://i.stack.imgur.com/5i3W0.png) What my boss wants is a way to input the population and the age and get out a predicted revenue number. What would be the best way to go about doing this? Should I combine every markets data together and then create an ARIMA model using population as an external regressor? Would that work? Or should I create separate models for every market? Is ARIMA even what I want to do? Maybe just a linear regression with both variables and test for interaction? I'm relatively new to professional statistics, having just graduated college in May, and any help would be greatly appreciated! edit: the dataset looks like this:               mkt    days.since.launch   date    revenue    mkt.population      for each market there are a few hundred datapoints.


title: GLM regression - help choosing model specification

text: I think I need to use a Poisson-family regression or negative binomial regression. My variables are as follows: Y is an integer value ranging from 0 to ~1200. It represents sums (number of species summed over an areal unit). There are in fact many zeroes but no negative values. X1 is a categorical variable, x2 is continuous (which also contains a few zeroes) and X3 is categorical. All are positive values. Variance of Y is larger than the mean.                        Y            X1         x2                X3        Min.   :   0.00   01:29551  Min.   : 0.000   2009 : 2474        1st Qu.:   5.00   02:72289  1st Qu.: 7.646   2010:28484        Median :  23.00             Median :13.000   2011:  882        Mean   :  77.21             Mean   :12.634                   3rd Qu.:  80.00             3rd Qu.:17.000                   Max.   :1155.00             Max.   :30.000                   ![Histogram of Y](http://i.stack.imgur.com/y1ala.png) Y is negatively skewed (i.e., skewed to the left). Histograms of residuals from a basic linear model (lm) and a QQ plot indicates the results are also skewed. The residuals plotted against fitted values also indicate that a linear model may not be appropriate because more points are above the line than below (across all values of x). Is it correct to use GLM with a poisson distribution with log link in this case?               Mydata.poisson  <- glm(Y~X1 +x2 + X3 +x2:X3, family=poisson, data=mydata)      Or more specifically, should I use the quasi-poisson? (in the regular poisson, my df was “31839 Total (i.e. Null); 31833 Residual”, Null Deviance was 1085000 and Residual deviance was 1079000). Also I believe this would be a case where I need to use a zero-inflated model? I am confused as to how to set this kind of model up. I read that a negative binomial distribution is similar to a poisson distribution, and better to use when the variance of your Y is greater than its mean, but isn't a binomial regression used when your response is binary? EDIT: I have used the following negative binomial model:               Mydata.nb  <- glm.nb(Y~X1 +x2 + X3 +x2:X3, data=mydata)      I understand that one should still check the residuals to see if the assumption of linearity holds (e.g., see discussion here: What are the assumptions of negative binomial regression?). A plot of the standardized residuals is included below and suggests that perhaps the relationship is not very linear. Would you agree? How can I resolve this? ![Residual plot for Negative binomial model](http://i.stack.imgur.com/dZDCM.png)


title: What is the proper terminology for the parts of a bayesian hierarchial model?

text: Lets say I have the following JAGS model, what are the proper terms for the three parts and 1/tau?                 model{       #Term???       for(n in 1:Ndata){       y[n] ~ dnorm( mu[ subj[n] ], tau[ subj[n]] )       }            #Term???       for(s in 1:Nsubj){       mu[s] ~ dnorm( muG, tauG)        tau[s] ~ dgamma( 5, 5)       }            #Term???       muG ~ dnorm( 80, .01)        tauG ~ dgamma( 1, 1)            }      What would I call each of the three parts? Also, for the parameters I would say:   * mu= "Individual Level Mean"    * tau= "Individual Level Precision"   * muG= "Group Level Mean"    * tauG= "Group Level Precision" Now what would I call 1/tau? "Variance"?


title: Calculating weighted Euclidean distance with given weights

text: This question is regarding the weighted Euclidean distance. I have three features and I am using it as three dimensions. I need to place 2 projects named A and B in this 3 dimensional space and measure the distance among them. But the case is I need to give them separate weights. if first feature is given weight 1, second feature should be given weight 2 and third feature must be given weight 3. If values of project A are (4,3,2) and values of project B are (6,2,7) can I use the weighted values as it is to calculate the distances between 2 projects? Please someone kindly help me in solving this? If this method is wrong is there any other way that I can use this? The basic thing what I need is giving them weights "what I need". This is just an example. I would be needed to expand this upto a n dimensional space. It would be able to be achieved using Euclidean distance. But main problem what I am facing is giving them weights what I am looking for.


title: Predictive analysis of data using R

text: I am trying to analyze data to predict the effect of 8 predictors on a response variable. I ant to do it using R. There are 8 predictor variables (x1, x2, ..., x8) and one response variable (y). I would like to find out which predictors have the largest influence on the response variable. I have already looked into models like logit regression and probit regression, but all of them need the response to be in 0/1 format. I am relatively new to R. Any help will be appreciated. Sample data:               x <- structure(list(S.no = 1:39, x1 = c(3.250596373, 2.450039885,      4.510736811, 2.429919206, 3.337809323, 1.492779301, 1.954807753,      5.154587813, 4.245647266, 2.826065366, 2.905970969, 2.60327762,      1.707142378, 1.082110181, 3.018482459, 2.036097259, 1.76077728,      3.390921555, 2.717923955, 4.526982605, 3.56285367, 2.067214422,      3.488982299, 3.562666419, 2.438279144, 4.238965876, 3.142895899,      2.594398665, 2.698489796, 3.309983984, 3.720543973, 3.67107841,      3.127048591, 1.363051041, 2.049245212, 2.777363784, 0.620092012767157,      3.51349138, 3.364191728), x2 = c(5.111023588, 4.969426778, 4.358447886,      3.898288664, 5.143436611, 4.82215687, 9.046313251, 4.632897736,      7.179105356, 4.789607667, 6.728836336, 3.819369146, 3.984776837,      6.26260718, 6.872900728, 10.19153143, 5.908181482, 6.811558747,      3.078856263, 6.522185906, 0, 4.048535963, 7.106468995, 8.295028612,      1.624650732, 4.183036858, 5.110667077, 2.377417101, 4.808793002,      1.975770073, 4.595705252, 4.341435526, 7.287675726, 5.305938153,      9.192668385, 7.919791579, 2.14634204, 4.279065244, 3.456415183     ), x3 = c(0.0908320519238175, 0.13266649971562, 0.258244321357115,      1.138276093, 0.0164448266821083, 1.464583483, 0.0960413836822101,      0.342061911817791, 0.119450732180513, 0.0290706457259754, 0.404738201124871,      0.203768610333403, 1.36915744, 0.0688836208234231, 1.211282801,      0.120083034038544, 0.314686886240134, 0.0156086673227241, 0.348031451601191,      0.598582980510635, 0.423262496464853, 0.0702598617576042, 0.198647673862676,      0.225821359083056, 1.824161863, 0.341704912159248, 0.0995473011086384,      0.190999270416796, 0.310999015642607, 0.485539790835786, 0.180986322152118,      0.77598437447802, 0.0429616530115406, 0.346864863595271, 0.434614085449358,      0.381001881620193, 0.333510559038226, 0.0250939546773831, 0.259975278336901     ), x4 = c(1.983757737, 20.50733268, 3.640023144, 3.528323863,      15.99619436, 7.489093737, 16.2062464, 9.816135606, 6.84983451,      2.912239466, 5.437503567, 35.8734624, 1.350936321, 28.55887031,      24.51023709, 49.29035933, 4.492149826, 122.0338801, 3.58950992,      17.11606415, 12.88815937, 5.312689839, 8.181119572, 36.36132304,      6.175496541, 5.078362324, 25.86563698, 8.01586315, 6.409855682,      8.085280012, 1.698834536, 3.809967237, 23.40496033, 8.041996294,      24.28479333, 2.051865304, 15.17903992, 10.56367934, 15.56299108     ), x5 = c(7.731494033, 5.673013527, 9.645267387, 12.47397681,      8.359512774, 4.737791104, 9.552434843, 4.67274363, 6.999587498,      6.463094738, 13.84710176, 14.15624423, 10.31230666, 11.28808186,      6.534183339, 12.40019882, 11.79614939, 12.38332647, 10.90812832,      7.321543633, 5.807470868, 9.854882389, 5.572464318, 10.28845247,      15.43592456, 15.30518284, 9.104886672, 14.25091879, 10.75910788,      7.648387336, 9.363347563, 4.304566946, 10.30984668, 8.647368931,      17.05992869, 10.34736023, 9.3230891, 7.322386746, 11.84261719     ), x6 = c(3L, 5L, 7L, 4L, 1L, 2L, 7L, 3L, 7L, 5L, 4L, 6L, 8L,      3L, 4L, 1L, 2L, 4L, 4L, 0L, 0L, 6L, 3L, 7L, 2L, 4L, 7L, 3L, 0L,      1L, 9L, 4L, 2L, 6L, 5L, 5L, 5L, 4L, 4L), x7 = c(3L, 2L, 4L, 1L,      4L, 5L, 3L, 1L, 1L, 2L, 3L, 3L, 0L, 4L, 5L, 1L, 4L, 4L, 6L, 5L,      3L, 2L, 5L, 3L, 0L, 2L, 3L, 2L, 5L, 1L, 1L, 2L, 5L, 2L, 3L, 4L,      0L, 7L, 0L), x8 = c(0.886405834295955, 1.055043892, 0, 1.019546273,      0, 0.363985182144308, 0.679987758960231, 0, 1.411360215, 0.990086072391923,      0.89255222931602, 0.0010880925135689, 2.303223251, 0.68924500650167,      0.861024867943333, 1.508676629, 0, 0, 0.633267197954753, 0, 0.352350421530178,      0.393923070064454, 0, 1.024464645, 0.577382403320939, 1.715032197,      0, 0.960435412665205, 0.806780641058154, 0, 1.466079364, 0.906711779133489,      1.965319112, 1.077518936, 0.538411084069564, 0.681301829171748,      0, 1.561002057, 0.124320537826079), y = c(12.40978817, 10.98435301,      12.82285819, 20.26438109, 12.95578413, 19.68995622, 31.19602354,      9.752948098, 19.01595709, 9.834666871, 23.69686387, 10.7471251,      19.41198687, 17.85098018, 28.66033204, 37.11149113, 19.34762966,      19.65906474, 12.95615703, 23.56218774, 5.458084618, 8.155950489,      22.94092426, 28.25918554, 29.03745489, 13.06245528, 12.61267288,      6.605952449, 16.05423534, 6.014818848, 9.73480676, 11.62296417,      23.14633388, 14.38379501, 38.58286932, 28.85006327, 4.644943375,      13.2872141, 7.061056339)), .Names = c("S.no", "x1", "x2", "x3",      "x4", "x5", "x6", "x7", "x8", "y"), class = "data.frame", row.names = c(NA,      -39L))


title: Equivalence test for two multivariate normal distributions?

text: I'm trying to compare two samples from multivariate normal distributions to see if their distributions are equivalent (within a factor of epsilon). The standard version of this test is the energy test, but it is not useful for my purposes because it uses $P = Q$ as the null hypothesis, whereas I need $P \neq Q$ to be the null. Some previous work has been done on this, mostly in this book, but nothing for my situation. How can I extend the interval inclusion method described in this book to use the energy statistic? Should the test statistic used to get the confidence interval be the original test statistic (i.e. based on the null hypothesis $P = Q$), or do I have to derive my own based on the new null hypothesis? I think derive my own is definitely beyond my capabilities at this point.


title: How to compare changes in group means over time

text: I would like to compare group means among 4 groups and analyse changes in their characteristics over time. I have a panel data set which is unbalanced. My sample sizes are not equal and therefore violate the assumptions of a mixed model anova? Any suggestions on how I can proceed? Also I wanted confirmation if a mixed model anova is the same as a multilevel modelling.


title: Structural shifts

text: I am using the R programming language to analyze the cointegration of 5 data sets. I have determined that the series share 3 cointegrating vectors after taking structural shits into account following the approach of Lutkepohl et al (2004). My questions are:   1. How can I incorporate the 3 vectors into the full matrix form of the final VECM?   2. How do I code in R the estimation of a VECM with a shift dummy variable to capture the break detected by the Lutkepohl et al (2004) test i.e., is it possible to include a shift dummy in a VECM? If so, can you please help me with a code example; e.g say the five data sets are `a b c d e`, cointegration rank `r=3` and break date `1975Q1`. Any references would also be very much appreciated.


title: Relative entropy optimization

text: I have what I hope should be a relatively quick and easy to answer question. Would it ever make sense to replace the probabilities by an objective function within the relative entropy formulation? I imagine a form of robust decision making under which a decision maker wishes to choose some value $s$ that maximises some function $g(x,s)$ say, but is uncertain as to the underlying distribution of $X$ except that is has distribution $F \in \mathcal{D}$ say, where $\mathcal{D}$ may be a constrained set of distribution functions (e.g. moment constraints). They will want to choose the value of $s$ that is most 'robust' against the distributional uncertainty relative to the true optimum, call it $g(x,t)$. This may look something like: $\max_{s\geq 0}\min_{F\in \mathcal{D}}\min_{t\geq 0}\int_{\Theta_X} f(x)\log\frac{g(x,s)}{g(x,t)}\,\mathrm{d}x $ In words the problem itself might be though of as the average of the logarithmic relative regret of choosing $s$ over $t$. If anyone has any thoughts on why this might be a bad idea, or problems that I would likely run in to then I would be grateful to hear from you. My main question comes really is whether this would make sense?


title: Linear regression optimization

text: I have the following linear model with 2 predictors $x_1, x_2$ and a cubic transform on $x_1$: $$ \hat{Y} = c +\beta_1x_1 + \beta_2x_1^2 + \beta_3x_1^3 + \beta_4x_2 $$ Where $x_1$ can assume negative values. I ran several regressions on subsets of the original population based on the absolute value of $x_1$ and noticed that $R^2$ improves gradually as I increase the minimum required abs value of $x_1$. For example, $R^2$ is higher for the sub-population where abs($x_1) > 10$; sub-population abs($x_1$) > 20 has higher $R^2$ then abs($x_1$) > 10 and so on. Can I add an interaction between abs($x_1$) and the full set of terms in the model? Or are there linear optimization techniques to solve this? **UPDATE** By imposing a constraint on the _absolute_ value of $x_1$, effectively we are saying the more the information available, or the stronger the signal, the higher the quality of the prediction. I wonder if this is still a case for interaction analysis, or if there are other techniques that can be used more effectively. Does it make sense to distinguish between the non-linear relationship between $x_1$ and $\hat{Y}$ represented by the cubic transformation (let's call it a "quality relationship"), and an "information quantity" term represented by the absolute value of $x_1$?


title: Set of numbers statistically higher than other

text: I have a set of numbers (number of publications for each year of one person). I don't know the exact distribution. Now I also have another set of numbers to compare that to (number of publications for each year of the same person during another time period). I'm interested whether the person has a higher publishing rate in that second time period. Can someone point me to a simple method I should use to determine the probability that one set of numbers is usually higher than the other? (ideally a source where I can find all the maths) Does it make sense to first sum up all publication numbers per year to get an overall for each period or do I lose information this way? Can I get more precise results if I assume a known distribution? Any suggestions which distribution makes sense for publication rates?


title: Subtracting t-distributions

text: I have two linear regression parameters of interest, b1 and b2. Both parameters are from a linear model built from 14 datapoints and having 7 model parameters, including an intercept. Interest lies in testing the null hypothesis b2 - b1 = 0 versus its two-sided alternative. I have two questions concerning this test: First, according to Slutsky's lemma, the subtraction of model parameters b1 and b2 converges to the subtraction of the distributions of b2 and b1 (t-distributions with n-p=7 degrees of freedom). However, I am afraid I am not allowed to work asymptotically since I have such a low number of data points, or can I? Second, how can I calculate a subtraction of these two t-distributions? Thank you, Koen


title: Comparing different types of losses as functions of lambda?

text: The usual pictures we see when dealing with different loss functions look similar to this: ![losses](http://learning.cis.upenn.edu/cis520_fall2009/uploads/Lectures/hinge- exp-logistic-bounds.png) Here we see y*f(x) on the x-axis with an error associated with it. Suppose I have a logistic regression where I am trying to use cross-validation to find the optimal lambda. I generate a list of lambdas and then build a model for each value of lambda and check what the loss function would give me at that specific value of lambda. I then find the lambda which minimizes the CV error, be it misclassification, hinge, or deviance. Now my question is what these plots are supposed to look like as functions of lambda? The CV-Error I get changes with respect to lambda (obviously) and each loss function is on its own scale. I'm not sure if I have made a mistake because I don't know what scales the loss functions should be on and what shapes they normally take when optimizing for different parameters. The strange thing is that all 3 error types tell me that a lambda value of 0 (no regularization) minimizes the loss function. I know through another R implementation (cv.glmnet) that this is not the case. Either I'm doing the cross-validation incorrectly or I wasn't sure how to define my loss functions properly. Any help regarding any of these issues would be much appreciated.


title: Which statistical method to use for finding systematic patterns in data

text: As part of a broader study I am analysing 30 websites that fall into 3 categories:   * Consumer (10 sites)   * Commercial (10 sites)   * Health (10 sites) The approach I used was a 'tick and flick' spreadsheet with 24 dichotomous variables that represent features of the website that are either absent or not (i.e. they receive a tick if they exhibit that particular feature). Here is an example of the data. The numbers represent how many websites from each category contain each particular feature (variable). I want to know which kind of statistical test would be used to find if there are any systematic patterns about which 'Category' of website tends to correlate with particular variables. For example, which websites tend to share power with users to edit/contribute web content (measured by variables 2,3,4,5,6,7,8,13,16,19,23,24)? I would rather use a more robust/rigorous statistical approach than simply counting up totals, or 'eye-balling' patterns in the data. Thank you in advance.


title: Simple linear regression - understanding given

text: The question is to fill out the missing numbers (A-L) of a simple linear regression model. I am having problems with converting and interpreting the given table in terms of variables. Would it be possible for someone to confirm and clarify things for me. The first table represents regression statistics ![](http://i.imgur.com/qeT0u9m.jpg) True model $$ Y_t = \beta_o + \beta_1X_t + \mu_t $$ Estimated model $$ \hat Y_t = \hat\beta_0 + \hat\beta_1x_t $$ This is what I am confused about   * Does the first standard error (12.8478) mean $\sum\hat\mu_t^2$ ?   * Does the standard error for the intercept in last table (14.6208) mean $\sum\mu_t^2$ ?   * Does 3.8508 equal $\hat\beta_1$ ?   * In order to calculate RSS (for J) I need $\sum \hat\mu_t^2$ does this confirm that my first two points are incorrect   * I know $G=\hat\beta_1^2\sum x_t^2$, how do I find $\sum x_t^2$ If I am wrong, would it be possible to know what those numbers mean in terms of variables


title: When would I choose Lasso over Elastic Net

text: What are the scenarios where Lasso is likely to perform better than Elastic Net (out of sample prediction)?


title: Does zero correlation between 2 differenced series implies no cointegration between original series?

text: The question is related to this one. In this question @mpiktas gives an answer on why checking correlation is not enough but the answer doesn't seem completely correct to me for the following reason: If 2 time-series are cointegrated, i.e. there is a linear relation between them $$y_t = a + b x_t + \varepsilon_t$$ with stationary $\varepsilon_t$ it _implies_ a linear relation between their differences $$\Delta y_t = b \Delta x_t + \Delta \varepsilon_t.$$ So if series are cointegrated, their differences should be correlated. And this means that if we don't see significant correlation between differences in our data, then there is no cointegration either. Is this correct or am I missing something? The question arises because I look for relations between hundreds of time- series (mainly non-stationary) and the way I do this is by considering correlations between their differenced counterparties. And I assume that if I don't see correlation between differenced series there is no cointegration either - it suffices to check only correlations.


title: Pairwise vs. total independence of discrete uniform random deviates

text: > Let $X$ be a discrete uniform random variable on the set $\\{000, 011, 101, > 110\\}$ of four binary integers, and let $X_{i}$ denote the ith digit of > $X$, for $i = 1, 2, 3$. Show that $X_{1}, X_{2}, X_{3}$ are independent > pairwise, but not totally independent. Can you generalize this example to > more than three random variables? Can anyone help me with this exercise?


title: How to perform a t-test to compare treatment and control groups when the experiment has been replicated three times?

text: I have data from a experiment that has been repeated 3 times under the same condition. Each time, 5 measurements of each of 2 groups (a control and a treatment) were taken. I want to run a t-test on my data, to infer if there's a signification difference between the control and the treated group.   * How should I use the test?    * Should I take the 15 measurements for each group and feed it to the algorithm (I'm using instat graphpad)?    * Should I average the 5 measurements of each repetition and feed the algorithm the mean and SD of each repetition?


title: What does the residuals effectively mean in this model?

text: I'm doing sme test using R and linear models. This is the code I'm using:               > a     [1] 0.10 0.05 0.02 0.04 0.09 0.05     > b     [1] 0.07 0.02 0.08 0.02 0.01 0.06     > c     [1] 0.03 0.03 0.01 0.05 0.06 0.09     >      > r.a     [1] 0.10 0.15 0.17 0.21 0.30 0.35     > r.b     [1] 0.07 0.09 0.17 0.19 0.20 0.26     > r.c     [1] 0.03 0.06 0.07 0.12 0.18 0.27     >      > mod = lm(r.a ~ r.b + r.c + 0)     >      > mod          Call:     lm(formula = r.a ~ r.b + r.c + 0)          Coefficients:        r.b     r.c       0.8764  0.5231            > mod$residuals               1           2           3           4           5           6       0.02295527  0.03973251 -0.01561396 -0.01929937  0.03054829 -0.01912000      >      > mean(mod$residuals)     [1] 0.006533791      The three serie above ( **a, b, c** ) are the returns of three stocks, and `r.a, r.b, r.c` are the cumulative returns `cumsum(a) cumsum(b) cumsum(c)` Now, as you can see I have used the `lm` function to calculate the linear model, at the end of the code there are the residuals of this linear model. I need to understand what those number effectively mean. In other words... I used the returns of these three serie, so are the residuals the difference of the percentage of return ? I try to exlpain in with an example: Take a look at the third value of each series(cumulative returns) Example:               r.a[3] - 0.17     r.b[3] - 0.17     r.c[3] - 0.07      And the third residual is: **-0.01561396** Does it mean that there is a variation of 1% ? (`0.01 = 1% 0.1 = 10% 1 = 100%`) And the second question is: What is the relathionship between the residuals and the mean of the residuals? Example: \- the third residual is: **-0.01561396** \- the mean is 0.006533791 Can I say THIS (third) residual is around 1% far from the mean? With this I mean...IF those three stocks(the returns of those three stocks) are moving and the spread goes to (~0.006533791) the total return of this movement is 1% ? Thank you!


title: Why are samples within a cluster less informative than randomly chosen ones from entire population?

text: Please give me mathematical explanation if possible. And also in the book Kothari 2004, it says: > There is also not as much information in ‘n’ observations within a cluster > as there happens to be in ‘n’ randomly drawn observations. Can you also give me mathematical explanation for this? My second question is why is simple random sampling preferable to cluster sampling? I'm interested in the mathematical proof of betterness of the former in terms of randomness.


title: Labeled multidimensional sequences

text: I'm trying to find similarities in 3-dimensional sequential data. The sequences are 3-uples $(t,r,d)$ each sequence is generated by one subject during a 3-6 months period:   * $t$ is a task identification ranging from 1 to 50,   * $r$ is the result obtained by the subject in the execution of the task (0 to 100),   * $d$ is the difficulty level of the task from 0 to 4. Every subject executes a different sequence of tasks : the tasks are executed in any order and the number of tasks executed by each subject is different (from 50 to 400). Some subjects execute tasks once a week, other twice, other 3 times, some subjects do not execute tasks for example during one week or two and then continue executing, all of them in a non specific order. Before execution period each subject is evaluated in 6 skills, and after execution period he is evaluated again in the same 6 skills to measure potential improvement, so that each sequence is labeled with a skills 6-tuple (that can be summarized in just one binary label (Improve, not improve). I have about 300 subjects and 50,000 executions. My problem could be similar to Amy McGobern's paper > McGovern, A, Rosendahl, DH, Brown, RA, and Droegemeier, KK (2011). > Identifying predictive multi-dimensional time series motifs: an application > to severe weather prediction, Data Min. Knowl. Disc., 22: 232–258. but their word size seems too short. I have up to 50 tasks (I could work with the 15 most frequent but lots of missings appear) and my sequences are also quite different in length from 50 to 400 tasks executions. Maybe functional clustering Peng's paper may help: > Peng, J, and Müller, H-G (2008). Distance-based clustering of sparsely > observed stochastic processes, with applications to online auctions. The > Annals of Applied Statistics, 2(3): 1056–1077. However, I'm not interested in the specific time point when a task is executed, I can assume that subjects execute tasks at fixed time points and analyse the problem without considering the sparcitity of the tasks executions.


title: R choosing the right classification approach for $ transaction volume categories

text: Our customers are Merchants and use our online payment service. Before they started using our service, they indicated how much $ transaction volume they will make per year. However this turned out to be only 50% true. We have categories:               A customer: high volume ($1.000k <)     B customer: medium ($500k - $1.000)     C customer: low ($500k >)      Turns out, that if Merchants claim they will make $ 1 Miollion a year, its only true 50%, usually they make less. So we have                Merchants claiming to be A but are B,C      Merchants claiming to be B,C but turn out to be A      Data like:                ID, transactions_$_year, trx_rating_pre, trx_rating_real, LOTS_OF_OTHER_DATA...      1, 1,500,000.00, B, A...      2, 150,000.00, A, C...      3, 100,000.00, A, C...      4, 50,000.00, A, C...      Now I want to create a classifier that learns from the data so that we will identify the A customers BEFORE they actually make $ transactions. We already get 50% right by asking. We have collected like 50 data variables (both quant and qual) from which we think they will help us predict A better than 50% like                Alexa Ranking when they signed up      Faceook likes      Twitter followers      Country      Legal form      Type and number payment methods used on website      Webdesign evaluation      Various other data regarding their website structure      So overall this is just a binary prediction modell I think: predicting A customers vs. non-A customers. What package in R would you use for a learning model? Which model would you use and why ? Data Mining Models WIKI                e1071      kernlab      klaR      svmpath      shogun


title: When to use bootstrap vs. bayesian technique?

text: I have a rather complicated decision analysis problem involving reliability testing and the logical approach (to me) seems to involve using MCMC to support a Bayesian analysis. However, it has been suggested that it would be more appropriate to use a bootstrapping approach. Could someone suggest a reference (or three) that might support the use of either technique over the other (even for particular situations)? FWIW, I have data from multiple, disparate sources and few/zero failure observations. I also have data at the subsystem and system levels. It seems a comparison like this should be available, but I've had not luck searching the usual suspects. Thanks in advance for any pointers.


title: Syntax for means and standard deviation on repeated measures GLM in spss?

text: This may be a silly question, but when running a **repeated measure** GLM in SPSS, and entering syntax for means tables (/EMMEANS), is there any way to get that table to show not just the standard error, but the standard deviations as well, for all cells and interactions? I need to report the standard deviations with means for a paper. Or maybe there is a better way to get them? /PRINT=DESCRIPTIVE only gives sd's for for the DV at each independent level but not the interaction between them (e.g. for a, and b, but not a*b). Any best way to do this in syntax? Other ways?


title: Using plm to account for non-individual fixed-effects in R

text: This is almost dumb and a little embarassing, but I can't figure out how to (or even whether it is possible) to use the plm package in R to run a regression including fixed-effects that do not correspond to the individual observation unit. For instance, I have observations on firms' outcomes, but I need to include sector-specific dummy variables. Of course I would have to "demean" these out (because there are so many of them). Is that even possible?


title: Formulating a hypothesis test about dependence

text: I have a data set of 91 variables. It is related to a qualitative analysis I made to analyse a few phenomena. For this reason, the questions are grouped in 10 clusters, since they express different aspects of the analysis. I made a subset using the following R command:               example1<- subset(data, select=c("a1", "a2", "a14", "a21"))      then I made a dependence matrix to check the dependence among the single variables in each cluster:               p1 <- stat1 <- diag(ncol(example1))          colnames(p1) <- rownames(p1) <- colnames(exemple1)          colnames(stat1) <- rownames(stat1) <- colnames(example1)      rn <- rownames(p1)      cn <- colnames(p1)          ###loop for the p-values     for(i in 1:ncol(example1)){          for(j in 1:ncol(example1)){              a <- example1[, rn[i]]              b <- example1[, cn[j]]              r <- chisq.test(a,b)$p.value              p1[i, j] <- r          }     }          ###loop for the statistic     for(i in 1:ncol(example1)){          for(j in 1:ncol(example1)){              a <- example1[, rn[i]]              b <- example1[, cn[j]]              r <- chisq.test(a,b)$statistic              stat1[i, j] <- r          }     }          ### placing the p-values in the upper diagonal of stat     stat1[upper.tri(stat1)] <- p1[upper.tri(p1)]      diag(stat1) <- 1      stat1  # this is for the first subset      But now I have two questions:   1. How can I make the dependence hypothesis test among the variables and set my alpha?    2. How can I make the regression analysis?


title: Can I include year as an independent variable in this nested logistic regression design?

text:
I want to make a nested logistic regression in R with the package mlogit. I would like to test how producer's decision to enter organisations (14 organisations) or not is affected by different factors. Producer can be in specific organisation and other according to year. I was advised to use year as separate variables in columns. This is OK. But for my different organizations will it be the same? So if my producer appear in my database several times (6 times) in row, isn't it redundant? Have you an Idea of the struture of database that I can adopt? Here is how my data looks like in R:                  Year   Organisation   Member AGE.Member ...     1  2005 Organisation 1 Member 1         37     2  2005 Organisation 1 Member 2         32     3  2005 Organisation 3 Member 3         32     4  2005 Organisation 4 Member 4         35     5  2005 Organisation 2 Member 5         33     6  2005 Organisation 3 Member 6         33          'data.frame':   18 obs. of  4 variables:      $ ANNEE       : int  2005 2005 2005 2005 2005 2005 2005 2005 2005 2006 2006 2006...
          $ Organisation: Factor w/ 4 levels "Organisation 1",..: 1 1 3 4 2 3 2 3 3 2 ...      $ Member      : Factor w/ 9 levels "Member 1","Member 2",..: 1 2 3 4 5 6 7 8 9 1 ...
          $ AGE.Member  : int  37 32 32 35 33 33 32 32 33 37 ...


title: ANOVA Theorem Explanation

text: I am trying to figure out why the following holds: Given $y_{i}=E[y_{i}|X_{i}]+\epsilon_{i}$ that $E[\epsilon^{2}_{i}] =E[E[\epsilon^{2}_{i}|X_{i}]] = E[V[y_{i}|X_{i}]]$ Specifically I am trying to understand why $E[\epsilon^{2}_{i}|X_{i}] = V[y_{i}|X_{i}]$? Clearly, I need a refresher on conditional variance and the rules for expected values....


title: Definition and Convergence of Iteratively Reweighted Least Squares

text: I've been using iteratively reweighted least squares (IRLS) to minimize functions of the following form, $J(m) = \sum_{i=1}^{N} \rho \left(\left| x_i - m \right|\right)$ where $N$ is the number of instances of $x_i \in \mathbb{R}$, $m \in \mathbb{R}$ is the robust estimate that I want, and $\rho$ is a suitable robust penalty function. Let's say it's convex (though not necessarily strictly) and differentiable for now. A good example of such a $\rho$ is the Huber loss function. What I've been doing is differentiating $J(m)$ with respect to $m$ (and manipulating) to obtain, $\frac{dJ}{dm}= \sum_{i=1}^{N} \frac{\rho'\left( \left|x_i-m\right|\right) }{\left|x_i-m\right|} \left( x_i-m \right) $ and iteratively solving this by setting it equal to 0 and fixing weights at iteration $k$ to $w_i(k) = \frac{\rho'\left( \left|x_i-m{(k)}\right|\right) }{\left|x_i-m{(k)}\right|}$ (note that the perceived singularity at $x_i=m{(k)}$ is really a removable singularity in all $\rho$'s I might care about). Then I obtain, $\sum_{i=1}^{N} w_i(k) \left( x_i-m{(k+1)} \right)=0$ and I solve to obtain, $m(k+1) = \frac{\sum_{i=1}^{N} w_i(k) x_i}{ \sum_{i=1}^{N} w_i(k)}$. I repeat this fixed point algorithm until "convergence". I will note that if you get to a fixed point, you are optimal, since your derivative is 0 and it's a convex function. I have two questions about this procedure:   1. Is this the standard IRLS algorithm? After reading several papers on the topic (and they were very scattered and vague about what IRLS is) this is the most consistent definition of the algorithm I can find. I can post the papers if people want, but I actually didn't want to bias anyone here. Of course, you can generalize this basic technique to many other types of problems involving vector $x_i$'s and arguments other than $\left|x_i-m{(k)}\right|$, providing the argument is a norm of an affine function of your parameters. Any help or insight would be great on this.   2. Convergence seems to work in practice, but I have a few concerns about it. I've yet to see a proof of it. After some simple Matlab simulations I see that one iteration of this is _not acontraction mapping_ (I generated two random instances of $m$ and computing $\frac{\left|m_1(k+1) - m_2(k+1)\right|}{\left|m_1(k)-m_2(k)\right|}$ and saw that this is occasionally greater than 1). Also the mapping defined by several consecutive iterations is not strictly a contraction mapping, but the probability of the Lipschitz constant being above 1 gets very low. So is there a notion of a _contraction mapping in probability_? What is the machinery I'd use to prove that this converges? Does it even converge? Any guidance at all is helpful. Edit: I like the paper on IRLS for sparse recovery/compressive sensing by Daubechies et al. 2008 "Iteratively Re-weighted Least Squares Minimization for Sparse Recovery" on the arXiv. But it seems to focus mostly on weights for nonconvex problems. My case is considerably simpler.


title: What's the convergence rate when solving L1 regularized optimization via coordinate descent with tiny step?

text: Wondering if there is an established result for the convergence rate when solving L1 regularized optimization via coordinate descent with tiny step? By "tiny step" I mean the step is always set to a very small positive constant, involving none of those step selection technique.


title: How to adjust for group size when analyzing interactions among group members?

text: Consider interactions between members of two distinct groups. If one group has substantially fewer members than the other, then interactions between that group and the larger group will be naturally more frequent than interactions between two members of that same group. Is there a way to take this into account when analyzing the interactions for frequency and duration of the interaction? Not quite sure where to start with this question. **EDIT 1** I have two distinct groups which can interact with each other. One of the groups is a lot larger than the other group. The members of the group can interact with each other or with members of the other group. I'd like to know if either group particularly favours interaction with the other group and if interactions between members of different groups last longer than interactions between members of the same group. As one group is a lot larger than the other group, if a certain member had an equal chance of interacting with ANY member of ANY group then I would observe a lot more interactions between members of the small group with members of the large group than interactions between members of the small group. This would possibly lead me to the conclusion that interactions between the groups are more likely to occur, even though we said at the start that there is a completely even chance of interactions occurring between any member of any group. How do I factor this into calculations? **Edit 2** The data is yet to be collected. I will be able to record the type and duration of any interactions as well as if the interaction occurred between groups or within either group. I can not identify individuals, only which group they belong to. Initially I'm interested in the frequency of interactions, but further down the line it might be useful to look at the duration of the interaction too. I plan to extend conclusions to a larger population. I plan to compare the results from these 2 groups to results I get from a different 2 groups. I'd like to try and identify factors influencing the willingness to interact between groups rather than within the group.


title: Numerical and graphical summaries of data distribution and relationships

text: I'm studying regressions. I was asked to write graphical and numerical summaries of the distribution of four variables and their relationships. What is a proper way to do it say for example in R? The data points are of the form $(x_i,y_i,z_i,t_i)$ where $0.5<x_i<35,\quad 0.3<y_i<14,\quad 0<z_i<9.15,\quad 0.4<t_i<19.$


title: Why don't we look at $R^2$ when fitting an autoregressive model?

text: $R^2$ measures explained variance. In an autoregressive model like `AR(k)`, we are carrying out a linear regression, and as such we would have an $R^2$ and an adjusted $R^2$. Why are they not used in practice?


title: Follow up of cluster analysis with membership prediction

text: I have 11 scale parameters for each of 218 observations belonging to subjects, I did standardized PCA to reduce dimensionality of the data and found two meaningful components. Using Euclidean distances this was followed by cluster analysis of these two components (explaining about 75% of the variance) with bottom-up approach using the hierarchical agglomerative clustering (HAC) by `FactoMineR` R package and Ward's linkage method. The optimal number of clusters was 4 as suggested by the package based on minimizing the ratio of two successive partition inter-clusters inertia gains.   This is just the number of observations per cluster:               > table(df$clust)            1   2   3   4        6  21  46 145      These 4 clusters turned out to be clinically important and subjects with cluster 1 were severely affected by disease. Cluster 4 were non-reactive subjects, Cluster 3 showed some reaction, and finally cluster 2 was like a special entity protected from disease. I don't know if these clusters can assume some kind of ordinal ranking or not. It is difficult to judge from the theoretical point of view related to the field, but I can say that cluster 4->3->1 is somehow showing some direction, and hence could be regarded as ordinal, on the other hand, cluster 2 is a little bit different but very important as subjects with this clusters were protected from disease. So, I am really confused as whether to consider these 4 clusters ordinal or not. Suppose that I have another set of 11 new readings of the scale parameters for one subject as new data, what statistical analysis would be useful to predict the membership of this subject to those 4 clusters? Could you please refer to a similar example with R code if possible? that would be greatly appreciated. Providing a professional answer would be highly esteemed, but also recommending some books using R code would also be encouraged, as I am searching for such a book that covers this topic thoroughly, many books are out there but it is difficult to judge which one would do the job. May be someone, has more experience with this kind of problems and can give a word of advise here.


title: Evaluate/Compare multiple numerical dataset for consensus or difference in R

text: I have a couple of stacked histograms which I need to compare/evaluate for similarity or difference. ![Stacked histogram generated from one dataset](http://i.stack.imgur.com/JK2GG.png) I believe rather than evaluating histograms is will be east to work with dataset used to plot these stacked histograms, which is in format:               RED                PURPLE       BLUE        GREY       YELLOW          22.0640569395   16.9483985765   0   60.987544484    0          8.1850533808    8.8523131673    0   82.962633452    0          6.8505338078    6.8950177936    0.756227758 85.4982206406   0.5338078292          6.7615658363    5.2491103203    1.6459074733    86.3434163701   0.6672597865          5.8274021352    7.384341637 2.1352313167    84.653024911    1.1565836299          7.8736654804    6.628113879 1.5569395018    83.9412811388   1.2010676157          7.1619217082    8.1850533808    1.2455516014    83.4074733096   1.3790035587          5.5604982206    10.2758007117   1.0676156584    83.0960854093   1.0231316726          7.1174377224    7.6067615658    0.7117437722    84.5640569395   0.756227758          7.8736654804    3.9590747331    0.6672597865    87.5    0.3113879004          7.6512455516    7.8736654804    0.5338078292    83.9412811388   0.5338078292          7.6067615658    8.9857651246    1.4679715302    81.9395017794   0.3558718861          8.9412811388    8.0071174377    1.3790035587    81.6725978648   0.5782918149          19.0836298932   9.2081850534    2.1352313167    69.5729537367   1.3790035587          14.9911032028   11.0765124555   3.2028469751    70.7295373665   1.0676156584          15.3914590747   10.8985765125   3.024911032 70.6850533808   1.2900355872          17.4822064057   12.5444839858   2.4911032028    67.4822064057   1.334519573          15.8362989324   13.0338078292   2.0017793594    69.128113879    1.334519573          17.037366548    10.4537366548   2.4021352313    70.1067615658   1.2010676157          20.2846975089   10.0088967972   0   69.706405694    1.0676156584          28.7366548043   12.6334519573   0   58.6298932384   0      Is there any possible way I can compare such dataset from multiple experiments (n=8) and visually show (plot) that these datasets are in consensus or differ from each other? Can QQplot be used to evaluate such datasets. Awaiting reply, Atul


title: Seeking a real illustration of Fisher's quote about DoE

text: My team and me would like to give a presentation to the non-statisticians of the company about the utility of the design of experiments. These non- statisticians are also our clients and they usually don't consult us before collecting their data. Do you know some real examples which would well illustrate Fisher's famous quote _"To call in the statistician after the experiment is done may be no more than asking him to perform a post-mortem examination: he only may be able to say what the experiment died of."_ ? Preferably we are looking for an illustration in an industrial/pharmaceutical/biological context. We think of an example of an inconclusive statistical analysis which could have been successful if it had been preliminary well designed, but maybe there are other possible illustrations.


title: Comparing the distribution (or properties of distributions) of large datasets for practically significant differences

text: I have several populations with a large n (100,000 - 1M+ samples), which i wish to compare. I know that virtually any statistical test i apply will reject the null hypothesis due to the extreme power of the statistical tests to resolve small differences. I am able to work around this when comparing means by setting my null hypothesis as $$ H0: |\mu1 - \mu2 | > \delta $$ $$ Ha: |\mu1 - \mu2 | <= \delta $$ Where the $\delta$ is the practically significant difference i care about. However i am lacking an intuition to apply a similar sort of transform for comparing variances (or standard deviations) of two populations as the hypothesis $\ H0: \sigma1=\sigma2$ is tested as a ratio. My initial thoughts was to subtract $\mu$ from each population (To give them a mean of zero) and apply a Chi Squared Goodness of fit test to the population histograms, however even with several hundred bins the n of the population observations still results in rejecting the null hypothesis for near identical distributions. I know i can game the power of the Chi2 test by:   1. Subsampling the data;or    2. Pretending n_observations is smaller then it really is.      * i.e Normalise the Expected and Observed histograms and multiply both by an arbitrary n, which is smaller then the actual observed sample size. The first option seems like a waste, besides this may be problematic as the data is obtained from many sites world wide. I'm not sure if the second is even valid, is anyone able to shed some light on this? How would i be able to justify my selected N? Is there a way to determine this by a desired practical effect size? Another alternative was to divide the variance of each of the observed populations by a control population, which at least allows me to rank effect size by % variance increase but then where do i go from there? I have ratios, and there is literature about estimating the error of a variance, however i will have no valid n to do a t-test of the difference of the ratios. Any suggestions welcome. * * * Sorry Whubar i accidentally hit enter and then took too long to finish getting back to you before editing closed, response below: We have a system taking measurements at multiple places within a site. Because of variation in local conditions, it is not possible to always use the same acquisition parameters for all measurements, some operators:   * have optimised the acquisition parameters per measurement location within a site;    * others optimise the acquisition parameters once per site;   * exclude the problematic measurements and base their data on nearby measurements; and    * have no needed to make any changes at all. We originally believed we could compensate for changes in acquisition parameters however when we blocked the data by acquisition parameters, we find the shape of distributions has changed marginally, which is expected, and a large mean difference. We don't really care about the absolute value between locations, just the delta between measurements at a particular site. We want to determine if we should:   * Encourage a independent location based optimisations and determine a fixed offset to subtract based on the optimisations applied; and   * Compare whether we are better off estimating missing values when the measurements are not available or subtracting a fixed offset based on optimisations;   * Specify that if optimisations are necessary that they go back and measure all new locations at the same acquisition parameters. Because we have lots of different sites, and we expect the data to be comparable between sites, we would expect the better method to have a lower variance, for a specific measurement location. However there are of course actual differences between the measurement sites, so we want to use the unoptimised measurement set as a control.


title: Distribution comparison across different treatments

text: if I have a table like this ![enter image description here](http://i.stack.imgur.com/ee6p8.png) In which, site = different sites of measurement (tissue or organ) type = different types of virus And I want to ask something about the data 1.What test can I use to test if the distribution of type is similar within site 1, site 2, site 3, etc. It looks similar, but I don't know to proof it statistically. To compare between two sites, I think we can use chi square goodness of fit test, but what if I want to do it simultaneously across all sites ? 2.Is there any ways to tell us that the distribution of type is dominated by only one type ? For example in this data set, we can observe that type 2 is the most dominant type across all sites. How can we do it statistically ? Thanks


title: How to assess normality in variables to determine whether to use Pearson's correlation or Spearman's r (or Kendall)?

text: I have one dependent variable and am trying to see if any of my 12 independent variables correlate with it, however I need to check everything for normality. I understand I need to look at skewness and kurtosis. Does the formula look like the following: skewness score subtracted by mean of condition divided by the standard error of skewness (likewise for kurtosis). Also are there absolute values that determine whether you should use a parametric/nonparametric test? E.g., skewness number above 2 = non parametric? I am using SPSS.


title: Probability that uniformly random points in a rectangle have Euclidean distance less than a given threshold

text: Assume we have $n$ points in a rectangular with bound $[0,a] \times [0,b]$, and these points are uniformly distributed in this plane. (I am not quite familiar with statistics, so I don't know the difference between uniformly pick a node in the area $[0,a] \times [0,b]$, or uniformly pick $x$-axis from $[0,a]$ and $y$-axis from $[0,b]$ independently). Given a distance threshold $d$, I may want to know the probability that two points' Euclidean distance is less than $d$, or more precisely, how many pair of nodes' distance will be less than $d$? * * * Maybe the following description would be unambiguous. Let me specify this problem. Given $n$ nodes and threshold $d$. These $n$ points are uniformly distributed in a rectangle $[0,a] \times [0,b]$. Denote a random variable $\xi$ as the number of pairs of points within distance $d$. Find $E[\xi]$.


title: Derivation for the confidence interval for a population proportion

text: I am having a conceptual hard time understanding where this formula came from. It does not seem to make any sense to me. Could someone shed some light on this: > The natural estimator of $p$ is $\hat p = \frac{X}{n}$, the same fraction of > success. Since $\hat p$ is just $X$ multiplied by a constant, $\hat p$ has > an approximately normal distribution. $E(\hat p) = p$ and $\sigma_{\hat p} = > \sqrt{\frac{p(1-p)}{n}}$. Standardizing, this implies that: > > $$P\left(-z_{\alpha/2} \lt \frac{\hat p- p}{\sqrt{p(1-p)/n}} \lt > z_{\alpha/2}\right) \approx 1 - \alpha $$ Could someone derive this equation so that it makes sense to someone who's never seen this before?


title: Expectation of an estimator?

text: When evaluating an estimator in a frequentist setting, using MSE and let say to compute the Bias of the estimator we compute the expectation of this estimator, are we supposing that the estimator has a probability distribution? doesn't this contradict the frequentist argument that a parameter is not a random variable?


title: How to adjust average rating for sample size on rating systems with more than two categories?

text: After reading How Not To Sort By Average Rating which deals with confidence interval for a Bernoulli parameter how would you extend it to more than two levels? For example: Items are scored between 1 and 5 (1 is worst 5 is best). What is the best way to adjust the average score per item in order to take into account the number of scores it received (one 5 score should not give it an average of 5!)? ![Shame on you Amazon!](http://i.stack.imgur.com/k7A0J.png)


title: How can I visualise overlapping time periods?

text: I want to draw a graph of sales over time. However, the time periods I have available to me overlap. For example:   * January 2012 to January 2013: 10 million    * April 2012 to April 2013: 12 million    * July 2012 to July 2013: 14 million    * October 2012 to October 2013: 12 million The obvious way to draw the graph is to do a simple line or bar chart, with the periods as labels on the y-axis: ![enter image description here](http://i.stack.imgur.com/3yZsn.png) (of course, the y-axis should begin at zero too!) But this doesn't feel quite right, because it doesn't really show what's really going on - it looks like a discrete set of points, not an overlapping set of periods. Is there a better way to visualise overlapping time periods?


title: Multiple Imputations exclude "not applicable"

text: I'm wondering how I should do the following in SPSS. I have a dataset with missing data (at random). Some values are blank, because the question was "not applicable" to that person (eg. questions about walking asked to people in a wheelchair). I assume I shouldnt want those values filled in with Multiple Imputations. How can I exclude those values? Thanks in advance


title: numerical measure for an ordered list

text: At the end of every semester I have to input my students' grades into a database. Typically this means I go one by one down my list of grades from my excel spreadsheet and manually fill the grade into the db application. One way to check I have done things correctly is to compare the mean and standard deviation of the grades from my excel sheet and the application. I realised that this would not capture the error of swapping over two grades: i.e. my excel sheet has the grades as: 4,4,4.5,6,4 and I put in 4,4.5,4,6,4. Is there an simple numerical measure that takes into account the order of the list such that a quick comparison of that gives me confidence that all is well? I suppose one could do some kind of Gödel numbering, but I see that as rather overkill.


title: Is there an alternative to paired t-test two multilevel ordinal variables? (similar to McNemar test)

text: Suppose I have two variables that measure some form of severity of illness on identical ordinal scale; one measured before treatment and the other measured after the treatment. I want to test if there is an improvement between treatments. Non-parametric correlations wouldn't be suitable; they would tell me if there is a relationship between the condition before and after the treatment.


title: Length of time series and likelihood estimation

text: In ARMA (with normal errors) model estimation, are there any empirical studies or tests to judge the minimum number of observations (length) of the time series that are required such that OLS is an acceptable estimation of the model (in contrast to an exact likelihood)?


title: Hold-out Validation vs K-Fold Validation?

text: To me, it seems that Hold-out validation is useless. That is, splitting the original dataset into two-parts (training and testing) and using the testing score as a generalization measure, is somewhat useless. K-Fold seems to give better approximations of generalization (as it trains and tests on every point). So, why would we use the standard Hold-out validation? Or even talk about it? Thanks


title: What are correlations among between test statistics in multiple comparison procedures

text: In the book Multiple Comparisons Using R Bretz et al write that their R package multcomp improves power above other multiple comparison procedures by "accounting for the structural correlations between the test statistics" (pg 41). The default option for this multiple comparison procedure implements this and is called "`single-step`."   My questions are 1)What do they mean by "structural correlations" (their explanation eludes me) 2)Is their a standard way to refer to or cite this approach. "Single-step" is what it is called in the call to their R function, but this is a general description of multiple comparisons. In my field people are mostly familiar with the usual Bonferroni corrections and I'd like to avoid the critical eye of reviewers and also point readers toward work that can help them understand what is going on.


title: What is the purpose of introducing nuisance parameters?

text: From Wikipedia: > a nuisance parameter is any parameter which is not of immediate interest but > **which must be accounted for in the analysis of those parameters which are > of interest**. Suppose originally the parameters are all those of interest. If later nuisance parameters are introduced, I was wondering what are some possible purposes of introducing nuisance parameters? Is introducing nuisance parameter to make the model identifiable? Thanks and regards!


title: Compare two Medians

text: I'm taking an intro to data visualization class, but I don't have a strong foundation in statistics. For now I'm trying to do something basic with bar graphs. I basically want to plot the median price of two types of products. I'm trying to see how much the condition of an item affects its price and comparing it with another type of product.               median     price     (normalized     in some way)                  |                       |                    |              A  B      |        B     [] []      |       []     [] []      |    A  []     [] []      |____[]_[]_____[]_[]______ condition (2 for now)            used      new      A is product type A. Let's say clothing for example. B is product type B. Let this one be guitar amps. A typical piece of clothing can sell from \$20 to \$300. A typical new amp can sell from \$100 to the thousands. As you can see, these products have wildly different price ranges. For some products, the price between used and new drops only a little, because wear and tear doesn't matter much. For others (like clothing or cars), the price drops a lot. If I just plot the actual median price, then only the comparison between used and new of the same product type would be useful. It would be hard to compare with the other category. What sort of calculations can I apply to the data in order to make this useful for comparisons? Am I approaching this in the wrong way? Data-wise, right now I have access to the individual prices of items with a specific selling condition. I'd appreciate any help.


title: How to get the diagonal elements of a covariance matrix from its sparse precision matrix?

text: I have a equation to solve $Ax = b$, where $A$ happens to be the precision matrix of a multivariate Gaussian distribution. I can use either direct solver or iterative solvers to get the $x$ vector. However, I also need to get the diagonal elements of the covariance matrix given by $A^{-1}$. I cannot directly invert the sparse precision matrix $A$ since the dimension is very high and its inverse will be full. It won't fit in memory. But I want to get the diagonal elements of the covariance matrix $A^{-1}$. How can I get it?


title: Are benchmarking and post-stratification equivalent?

text: Does "benchmarking" and "post-stratification" describe exactly the same process? The Wikipedia article on Statistical Benchmarking claims it is; however, the article on Sampling also mentions "poststratification" but doesn't link to "benchmarking". Perhaps there are other applications of benchmarking besides post- stratification?


title: Probability Plot (P-P plot): to what corresponds the straight line?

text: I want to show that my data cannot fit a given distribution (in my case it's the Pareto type II distribution). I did a P-P plot, i.e on the y-axis I plot my sorted experimental data and on the x-axis the inverse of the theoretical cdf. I added my Matlab code below in case my explanation is not clear. Most of the time, on the P-P plots, I see a straight line and it's said that if the point follow this line you can assume that your data are distributed according to this distribution. But how do you plot this straight line? Is it just y = x? If so, why? If not, what exactly do we have to plot?               %%% Definition of the probability density function, the cumulative density function      %%% and the inverse density function     pdf = @(v,x) x(2)/x(1) * (1 + v./x(1)).^(-x(2) - 1);     cdf = @(v,x) 1- (1+v./x(1)).^-x(2);     invcdf = @(y,x) x(1).*((1-y).^-(1./x(2))-1);          %%% Function to minimize to find the parameters %%%     fobserved = @(x) -mean(log((pdf(islets_vol,x)))); % islets_vol = my experimental data     options = optimset('Algorithm', 'interior-point',...                        'MaxIter', 1000, ...                        'MaxFunEvals', 1000);          [xhat_obs,~] = fmincon(fobserved,[0.5;0.5],-eye(2),[0;0],[],[],[],[],[],options);                                            %%% P-P plot %%%     yvals = sort(islets_vol,'ascend'); % islets_vol = my experimental data     xvals = invcdf((1:numel(islets_vol))/(numel(islets_vol)),xhat_obs);      scatter(xvals,yvals);


title: Modeling a TAR model that can handle the problem of outliers in nonlinear data

text: Can somebody please share idea on how I can model a TAR model that can handle outliers in nonlinear data? I need to compare such model with the general form of TAR model. Which computer software can be applied? Thank you.


title: Bayesion priors in ridge regression with scikit learn's linear model

text: I'm using scikit learn's linear model to do ridge regression. Ridge regression penalizes parameters for moving away from zero. I want to penalize for moving away from a certain prior, with each parameter having a different prior. Is this possible with scikit learn's linear model? I know there's a BayesianRidge module there, but I'm not sure what it does.


title: STEP BY STEP approach to learn LDA (Latent Dirchlet Allocation)

text: I have an urgent need to understand Latent Dirchlet Allocation (LDA) for Topic Modeling. I tried several sources, but it seems I do not have required knowledge to learn this. My statistics knowledge is very limited. I know basic probability, conditional probability and a bit of probability distribution. Can any one suggest a STEP BY STEP approach so that I can understand LDA clearly? If you can break the steps into smaller parts ( can be primitive or very silly ), I will be happier. For example, you can tell me, learn Multinomial Distribution, or so on like that. Thanks for your help.


title: Error Propagation vs Standard Deviation

text: Say you have a set of data of lengths, widths, and heights of a rectangular box. Why do we have to use error propagation to calculate the std. dev. of the box's volume? Why couldn't we calculate the volume associated with each l/w/h datum, and take the std. dev. of that?


title: Are there frequent itemsets mining algorithms based on linear algebra?

text: Matrix equations of some sort or matrix decompositions often lie in foundation of many algorithms not only in data analysis but in many other fields. The wealth of already written libraries makes it tempting to formulate frequent itemsets mining problems in terms of said matrix equations and/or decompositions, yet after some search I haven't found any algorithms which would rely on linear algebra. Could some please show me if such algorithms exist or explain their absence?


title: Estimate of variance with the lowest mean square error

text: Regarding estimators of variance from a iid sample of size $n$, Karl Ove Hufthammer says Estimates of variance from an iid sample: > if they do have a normal distribution, dividing by n+1 (sic!) will give the > lowest mean square error. Why does it give the lowest mean square error? What does he mean by "(sic!)"? Thanks!


title: Does confirmatory = inferential, and exploratory = descriptive analysis / statistics?

text: Reading this webpage, I wonder:   * Are confirmatory analysis / statistics and inferential analysis / statistics the same concept?   * Are exploratory analysis / statistics and descriptive analysis / statistics the same concept?


title: How to prove that the chance of occurrence of event A increases in presence of event B, when both have a different Poisson distribution?

text: If we have the events:   * A is the number of rainfalls/month in a given period of time (rate of which is equal in each month)   * B is the number of hailstorms/month in a given period of time (rate of which increases in winters) Can we statistically prove the probability of event A increases whenever the probability of event B increases? If yes, how? EDIT: thank you for the valuable input. The reason I am interested in this question is because in EEG tracings for diagnosing epilepsy we look at "spikes" which are events occurring at more or less a constant rate (though the rate differs in awake and sleep, but I am only interested in spikes in sleep). Normally, the rate remains constant during sleep, but in few cases, it has been seen to increase in the presence of an event called "k-complex". It is hypothesized that these k-complexes actually increase the epileptogenisity and hence the "spikes". So if we have a poisson of "spikes during sleep" and "k-complexes" (rarer events during a specific sleep stage) during sleep, can we statistically prove that spikes increase when number of k-complexes are increased?


title: Creating statistically similar sequences

text: I have to run a test on a program, with some statistical tools to make sure its results are statistically acceptable. This program will take sequences of data information (a vector of integers or symbols with different length) and use some distance metrics (like Levenshtein distance) to find the distance matrix between these sequences and then cluster them. Now I want to simulate some information to see statistically, how the program is doing. What I have in mind is to create some sequences that are statistically similar and see how does the program cluster them. What option do I have here to do so?


title: Geometric distribution without replacement

text: On an attempt to solve this problem I've managed to reduce it to finding the expected number of white balls picked until one black ball is observed (let's call that value $v$). Except that, unlike the geometric distribution, this needs to be done without replacement. Hypergeometric distribution doesn't come to the rescue as the number of black balls picked is immaterial (and of course the white balls must be picked consecutively). So I am looking for a way to compute $v$, or at least its underlying distribution. I've tried to search online but couldn't find a distribution that readily fits this. So the problem is: > Consider an urn having total of $n$ balls, $w$ of them are white, and the > rest is black. Pick $j$ balls **without** replacement. What is the expected > number of white balls that will be picked before a black ball is observed ? > Notice that the total number of balls picked, $j$, is fixed. My guess is that cumulative hypergeometric distribution could be used, while dividing somewhere by the number of ways the black balls intersperse the white ones (after the first black ball is observed). But I have little intuition how to go through with this. I.e. to divide before taking the cumulative distribution or after? **Update** Based on a question in the comments below, to be clear: the desired value is the expected number of white balls observed until **either** a black ball is observed, or $j$ picks happened, whichever comes earlier. Thank you.


title: Is it a problem to have some uncorrelated factors in a CFA?

text: I have a CFA with four factors - they are all intercorrelated except for two. Is this a problem for the model to have two factors uncorrelated?


title: Cramer-Rao type bound for Information Gain

text: I am interested in the Bayes risk of some distribution $\pi$ $$ r(\pi) = \mathbb{E}_{\pi(x)}[ \mathbb{E}_{\Pr(y|d,x)}[L(x,\hat x(y|d))]], $$ where $L$ is some loss function and $\hat x$ is the posterior mean (here: Bayes estimator). Here $d$ is some experimental design parameter. If $L(x,\hat x) = |x-\hat x|^2$ (squared error), then the **Posterior Cramer- Rao bound** gives a lower bound on the Bayes risk. The Bayes risk is also identical to the expected posterior variance. So to optimally design experiments, I can use **negative expected variance** as a utility function to maximize. Now I would like to reverse this logic to find the effective loss function and a lower bound for **information gain** as an alternative utility function. That is, suppose the utility function I use is information gain. Then,   1. What is the loss function which this is optimizing?    2. What is a lower bound on the Bayes risk for this loss function?


title: Does a stationary process implies a normal distribution of the data?

text: My understanding is 'no', a stationary process does not imply a normal distribution of the data. However I haven't found a clear indication in my library or online. I am interested in other resources comparing stationary and heteroscedastic processes to data distributions. I understand a stationary process joint probability distribution does not change over time. But does it necessarily mean that the _shape_ of the distribution is normal?


title: How are error terms calculated for moving average model in R

text: For an ARIMA (0,0,1) model, I understand that R follows the equation: xt = mu + e(t) + theta*e(t-1) (Please correct me if I am wrong) I assume e(t-1) is same as the residual of the last observation. But how is e(t) calculated? For example, here are the first four observations in a sample data: 526 658 624 611 These are the parameters Arima(0,0,1) model gave: intercept = 246.1848 ma1 = 0.9893 And the first value that R fit using the model is: 327.0773 How do I get the second value? I used: 246.1848 + (0.9893*(526-327.0773)) = 442.979 But the 2nd fitted value given by R is : 434.7928 I assume the difference is because of the e(t) term. But I do not know how to calculate the e(t) term.


title: Minimum of two exponential variates: What's wrong with this derivation?

text: Suppose we have $\newcommand{\E}{\mathrm{Exp}} X \sim \E(\lambda)$, $Y \sim \E(\mu)$, and $W = \min(X,Y)$. I know that $W \sim \E(\lambda+\mu)$. I know how to derive it. But, I tried this alternate derivation that gave me a different distribution for $W$, and I still can't figure out what's wrong with it. I started with $\newcommand{\rd}{\,\mathrm d}\renewcommand{\Pr}{\mathbb P}f_W(t) = f_X(t) \Pr(X<Y) + f_Y(t) \Pr(Y<X)$ Now I need $\Pr(X < Y)$. Seems straightforward: $\begin{align} \Pr(X < Y) &= \int_0^\infty [1-F_X(t)] f_Y(t) \rd t \\\ &= \int_0^\infty e^{-\lambda t} \mu e^{-\mu t} \rd t \\\ &= \frac{\mu}{-(\lambda + \mu)} \left[ e^{-(\lambda + \mu)t} \right]^{\infty}_0\\\ &= \frac{\mu}{\mu+\lambda} \>. \end{align}$ And, if I do the same thing for $\Pr(Y < X)$, I get $\Pr(Y < X) = \frac{\lambda}{\mu + \lambda}$. So $\Pr(X < Y)$ and $\Pr(Y < X)$ sum up to 1 as expected. Encouraging. And now I substitute that into my original equation: $ \begin{align} f_W(t) &= \frac{\mu}{\mu + \lambda} \lambda e^{- \lambda t} + \frac{\lambda}{\mu + \lambda} \mu e^{- \mu t} \\\ &= \frac{\lambda\mu}{\lambda + \mu}\left(e^{-\lambda t} + e^{-\mu t}\right) \>. \end{align}$ That... That's no exponential. Where did I go wrong?


title: Using Conditional Optimizations in Training Algorithm of Linear Support Vector Classifier

text: I am trying to come up with an easier algorithm than usual nonlinear optimization algorithms for linear support vector classifier. My idea is   1. Once classifier is determined, it's easy to decide support vectors (data misclassified + data near the discriminant plane)    2. Provided that support vectors are fixed, it's straightforward to compute linear classifier based only on those support vectors. So my algorithm to compute Linear Support Vector Classifier is repeating following steps until convergence.   1. Based on given support vectors, train a linear classifier    2. With linear classifier computed in #1, find new support vectors. Due to the fact this is a pretty obvious idea, I think somebody must've done something like this already. I am interested in how to prove that algorithm actually lead to convergence. Does anybody have any reference or clue?


title: Assuming the true mean

text: I've collected fluorescence data from some bacterial cells. Each cell has a gene in it which can be induced to fluoresce. However, even without being induced, the gene will still fluoresce a little. I have subjected the cells to a series of concentrations of inducer, including 0% inducer, and measured the resultant fluorescence over time. The fluorescence produced from cells that have not been induced is, effectively, 'background' fluorescence. Therefore, I would like to subtract the background fluorescence from the fluorescence measured in cells that have been induced in order to determine the fluorescence that is caused by induction. $n = 3$ for all measurements. My question is, can I take the mean fluorescence of the uninduced cells and assume that it is the true mean? I.e. assume that its $SE = 0$. If I can do that, when I subsequently subtract the mean background fluorescence from the mean induced fluorescence, the final calculated $SE$ will be lower. Here's some example data:                Inducer:       0%              5%      1             3411            4965      2             3427            4784      3             3412            5379      Mean          3416.666667     5042.666667      SD            7.318166133     249.0385958      Given that the background will change for each experiment, I assume that the mean background measured is the true mean for that experiment. Is that not correct? I must state that my maths/stats skills are rudimentary at best so please bear that in mind when answering.


title: Standardization of data during multilevel analysis

text: Recently I obsessed about standardization of scores. I am looking into statistical literature to see if standardization of data is OK or not. a brief about standardization: http://en.wikipedia.org/wiki/Standardizing by doing the standardization we shift the mean to ZERO ( 0 ), and then the data scores are the number of standard deviation an observation is above or below the mean. My question: I have done a repeated measure study, over 32 participants and 20 observation each participant. dependent variable is electrocardiography (EMG). specificity I am interested about the standardization of data per individual. I can just analyse the raw data I already have without any manipulation, or I can first standardize data per individual and then do any analysis on it. people who usually do this type of standardization tell they are not interested to between-people variation and because EMG is a very variable data type, so this is better to standardize data before. I am still not convinced with this rationality. As I understood multilevel analysis can control for Between-subject variability and within-subject Between-observation variability separately. In this case I am hesitating about standardizing my data before analysis. to be more clear: I have one dependent variable (response) and one explanatory variable.


title: STATISTIC - COMBINATIONS

text: Well, I trying hard with this problem...but two of my answers are wrong. Please check it. A painter has 4 benches to paint. He has paint of the following colours: RED, YELLOW, WHITE and BLACK. In how many ways can the benches be painted if (a) he chooses the same colour Answer 4 (b) at least one bench is of different colour I've got 15 but the answer says 31??????? (c) all the benches are of different colours Answer 1 (d) at least 2 benches should be painted red I've got 11 but the answer says 10????????? THANKs REPLY as soon as possible please Ok bye


title: Which post-hoc comparison to use and why

text: **Hopefully this is now a much clearer picture of the question I want to answer:** For each of 60 different plans I ran 35 simulation replications of my Agent Based Model (ABM) simulation. Then for each simulation replication of each plan, I generated a score ($\geq0$), so that I ended up with 35 scores for each of the 60 plans. I now want to find the best plan which I define as the plan with the smallest mean score. Of course I can compute the mean score for each one of the 60 plans, and I can rank each plan according to the mean score. At this point how do I test whether the 1st ranked mean score is statistically significantly different to the 2nd ranked score. And if I find it is not, how do I test whether the first two ranked mean scores are statistically significantly different to the 3rd ranked score, and so forth. Eventually I would like to end up with a set of plans (although the set can contain only one plan) for which I can state that the set of plans are equally best (as far as being statistically significantly the same). **This is what I was considering doing:** For each plan I compute the mean score and then select the 5 plans with the 5 smallest scores. I then use the Kruskal-Wallis H-Test to determine if any of the 5 means are significantly different (since my samples do not satisfy the normality or homoscedasticity assumptions of ANOVA). If the p-value returned by the test is less than my significance level $\alpha$ then I know that at least one of the means is different and want to test for which one this is (are). At this point I need to use a post-hoc comparison test but am unsure as to which one I should use and why.


title: What's the relation between deep learning and extreme learning machine?

text: Often I have found deep learning and extreme learning machine discussed together. Based on my little knowledge of the subject my impression is that they are different methods with different aims. So what's the relation between deep learning and extreme learning machine ? Is there any important paper that I am missing and that I should read?


title: Find weight of features for feature selection

text: I have a data set of videos from which I need to recognize the emotion of the speaker. For that reason I have some markers on the face of the speaker. I detect their movement as the speaker speaks and for every frame find the change in the marker coordinates. There are a total of 65 markers (blue dots) on the face of a speaker so at the end of one frame I have 130 (x and y coordinates) stored. I have chosen the mean and standard deviation of those points over the video sequence as my features (based on a paper published in 2009) My question is from these 260 features per video (mean and standard deviation of each x and y coordinate) how do I reduce this to something more understandable or which produces a meaningful probability output(in my case an emotion). In the paper it used the plus l and take away r algorithm based on the Bhattacharyya Distance criterion but I just couldn't understand how to actually say one feature is better than another based on that criterion.