title: An estimator for a histogram/summation based calculation

text: I'm looking for some ideas on how I can estimate the results of a particular calculation on a data set. The data is a very simple mapping of X=>Y where X is a real value and Y is a positive integer. The calculation that needs to be performed is given value Q what is the sum of all Y having X < Q. What I want is a way to come up with an estimator of the sum for any given value Q. I'm looking for ideas, or pointers, to such an estimator. Some key points, applicable to the software in which it will be used:   1. The estimator must need a relatively small amount of data, preferably alterable such that is the number of estimator variables increase the quality increases    2. There are sensitive areas in the domain, for example the accuracy from A..B is more important than B+. In general the accuracy requirement goes down as we go further right in the domain.   3. Perhaps the most vital requirement, it must be possible to modify the estimator's data directly from the addition/removal of a data item in the original set (that is, I cannot calculate over the entire data set each time to get the estimator values -- though occasionally is okay) _(Update: If you have an idea but it doesn't meet this point, or aren't sure, please tell me anyway, I'm sure I could adapt many approaches to be valid for this point)_ My first approach is to treat this as a histogram and simply divide the data into buckets of fixed sizes. For those sensitive areas I simply have more divisions. Unfortunately I'm finding that this is requiring too much data for the estimator to be accurate on the lowest part of the domain. To improve the histogram approach I've been looking for curve approximations within the histogram buckets. Something like variance seems appealing since it can be done with Point 3. _This is kind of mix of topics together, from programming to statistics to algorithm design. I'm asking for ideas and direction here since I think the problem is most related to the data set / statistical aspect._ * * * Update for first comment:   1. A..B / B+ in Point 2 just represent an arbitrary part of the domain. For example between 0 and 5 accuracy is more important than in 5+. These values will be constants.   2. The data set is in constant flux via additions/removals. The data set at Time 0 may have an empty intersection with the data set at Time X. It is therefore important that modifying the estimator through a series of additions/removals results in the same values as if it were completely recalculated. Note when I say cannot be recalculated it would be okay if there is some way to do partial recalculation (think constant time in terms of algorithm complexity). Also to the add/subtract, in the areas of high accuracy it is assumed the addition/removal of a single element will be highly significant. Note that I'm also okay if you have an idea that doesn't qualify for this aspect, since I'm sure I can find someway to adapt it.


title: Calculating a smooth 90% limit for differences in a time series

text: I have 50 data sets. Each set has three related time series: fast, medium, slow. My end purpose is simple, I want to generate a number that indicates a relative degree of change of the time series at each point. That relative degree of change should range between 0-1 for all the time series and all the data sets. The scales of the data set range from .0001 to 100. To accomplish this, I calculate the differences in a time series, delta(ts)=ts(t) – ts(t-1). Now I am trying to calculate an upper limit of 90% of those deltas. In other words, draw a smooth line on those differences such that only about 10% of the differences exceed that line. I will use that 90% limit to establish a maximum to normalize the differences between 0-1. Is this the best way to do this? I’ve been working on this for months, mostly linear programmatic methods, with no success. And trying to get it to work across 50 data sets is killing me. I’m sure there has to be an elegant mathematical way to do this. I can’t be the first guy in town trying to normalize a relative degree of change of a time series. Any help or directions for research are greatly appreciated! Obviously my math skills are weak so examples would be most helpful. Thank you everyone for your time and brain power! (Not sure how to attach a spreadsheet so only include an image)![my time series and differences](http://i.stack.imgur.com/tKL4l.png)


title: What transformations preserve the von Mises distribution?

text: The von Mises distribution is entirely defined on the circle with a density given by $$f(x) = (2\,\pi\, I_0(\kappa))^{-1} \exp(\kappa \cos(x-\mu))\ ,$$ where $x$ is in an arbitrary real interval of length $2\pi$, $I_0$ is the Bessel-I function of order 0, $\mu$ is a location parameter, and $\kappa>0$ is a concentration parameter. My question is: What transformations of a von Mises random variable preserve the distribution family? A slightly simplified version of this question (considering an important special case) can be restated as follows: Is there a bijective function (defined on the circle) $T$ and a $\kappa>0$ solving $$(I_0(1))^{-1}\, \left|T'(x)\right|\, \exp(\cos(T(x))) = (I_0(\kappa))^{-1} \exp(\kappa \cdot \cos(x))\ ,$$ The answer is easy for $\kappa=1$, because $T$ might be chosen as the identity, a simple shift, or a suitable reflection, Thus, I am interested if this can be solved when $\kappa\neq 1$. Unfortunately, I am not aware whether there are any results available on this issue. This is a crosspost from Mathoverflow.


title: How to interpret results from R anova in quantile regression?

text: Here is a sample output:               anova(fit1,fit2);     Quantile Regression Analysis of Deviance Table          Model: op ~ inp1 + inp2 + inp3 + inp4 + inp5 + inp6 + inp7 + inp8 + inp9     Joint Test of Equality of Slopes: tau in {  0.15 0.3  }            Df Resid Df F value Pr(>F)     1  9     1337  0.5256 0.8568          Warning messages:     1: In summary.rq(x, se = "nid", covariance = TRUE) : 93 non-positive fis     2: In summary.rq(x, se = "nid", covariance = TRUE) : 138 non-positive fis      How to interpret the above results?? Does the `anova()` function give the best model, for `tau=0.15` vs. `tau=0.3`?


title: What's the best approach to statistically analyze announcements on the price of a stock?

text: I have stock price data from Yahoo and I have corporate announcements from various RSS feeds that give me the date of the announcement, the date of the quarterly statement, etc. What is the best approach I should take to statistically analyze how much, if any, a new announcement has effected the stock price, whether an announcement of a quarterly statement date ahead of time amplifies the changes, if any, that a quarterly statement would make, etc.?


title: How to bin a quantitative covariate for multiple regression?

text: I do a multiple regression in which spacing between lines is one of four variables. However, this spacing between lines varies so I have grouped them into five classes designated with numbers 1 through 5. The number 1 represents the case of no lines present; the number 2 represents that the spacing between lines exceeds 20 m; 3 represents a spacing from 16 - 20 m; 4 represents spacing from 11 - 15 m; and 5 represents a spacing of 10 m or less. Can the numbers 1 - 5 be called "dummy variables" or do they have to be called something different?. As I understand it, dummy variables can only have the values 1 and 0.


title: .999 signifiance query

text: I am getting the following result for one of my variables after running a logistic regression. ![enter image description here](http://i.stack.imgur.com/ofxSF.png) Basically I have a set of binary dependent and one independent variable. One of the dependent variables is never 0 when the independent variable is 1. I believe that is why I am getting the above result. How would you interpret this? The .999 indicates that it's not significant but surely it should be? Also, what would you say about the odds ratio? EDITED: Yes, my bad. I have one DV and about 10 IVs. All of the other IVs display normal ORs and 4 have significant p values. The 'endoc' one I included in the screenshot above is the only 'outlier'. Here is the rest of the regression table: ![reg](https://www.evernote.com/shard/s161/sh/452f8fe2-560e-479f-8570-046aca48dc7d/312c12ac09662b8f5ff3d13fe3719788/deep/0/Shaphiro_Score_Validation_plusdescriptivestats.spv --Document3----IBM-SPSS-Statistics-Viewer.png). and this is part of my data ![data](https://www.evernote.com/shard/s161/sh/85a4bb71-a61c-4500-9edd-7708268911b6/2e42154bf3a5104e3fba30452b9c5d82/deep/0/Shaphiro_Score_Validation_2.sav --DataSet1----IBM-SPSS-Statistics-Data-Editor.png) I am basically writing a report for this, so I just need to know what to say about this occurence?


title: Bayesian approach to testing an unfair die

text: Suppose you have a six-sided die that you suspect is not fair and toss it N times. What would be a Bayesian approach to estimating the probability of the six outcomes given that you suspect the die may be unfair? I know that I have to come up with a prior density and a likelihood in order to define what the posterior would be. I was thinking of using the multinomial as the likelihood, but I am not sure what to use for a prior since the die is suspected to be not fair. Any help is appreciated... Thank you!


title: Dummy variables in a multiple regression

text: I have a linear regression model with 3 independent variables (let's say `A1`, `A2`, `A3`) and 2 different dummy variables, one for the gender (`d1`) and the other one for the location (`d2`). When I estimate the model with all the variables included, some of independent variables are not significant, but when I add just one of the dummy variables, all of the independent variables are significant. Could you please let me know why? And when I estimate the model in the form of f(`A1`, `A2`, `A3`, `d1`), I get different coefficients for the independent variables in comparison with the ones for f(`A1`, `A2`, `A3`, `d2`). Why? In general, should I add all the possible dummy variables to the model at once, or it is possible to have different similar models with just one of the dummy variables included? Thanks.


title: Point-biserial correlation for non normally distributed data

text: Can I use a point-biserial correlation if my data is not normally distributed and my sample size $N = 120$?


title: extremely left-skewed response variable - how do I model this dataset?

text: This is a histogram showing my response variable. ![histogram of response variable](http://i.stack.imgur.com/yJZaQ.png) The response is # (or proportion? or percent?) of aphids eaten off of cards in fields, to model predation by natural enemies. Predictors: fixed effects are both categorical (i.e. crop type, season) and continuous (i.e., landscape variables - amount of land in production, mean field size), and the variable 'landscape' is a random effect. I've been working on this so far with the response variable reflected ('number of aphids remaining' rather than 'number of aphids eaten') so that it's right skewed instead which seemed more possible. But I would rather work with it as left-skewed if possible, the results will be easier to discuss that way. Transformations do not help the response variable be less skewed. GLMM with Poisson errors does not work because the models created are too overdispersed. GLMM with negative binomial errors - same problem. Any way I try to model this data that treats the response variable as continuous, gives me problematic things when I'm model checking, like residuals by fitted plots that have clear patterns in them, and overdispersion. I'm starting to wonder if I need to somehow rank the data in the response instead? Or partition it into categories? The potential solution I've gotten to is hurdle models: "Hurdle models partition the model into two parts: a binary process generating positive counts vs. zero counts, and a process generating only positive counts. The binary process is modeled using a generalized logistic regression, and the positive count process is modeled using a zero-truncated count model" (paraphrased from Zeilis, Kleiber & Jackman 2008) Is there a way to do this in R? Or would I just do the two models separately, and discuss them separately? Or is there a way to get an AIC value for a hurdle model? Does anybody have any other ideas for how to model this dataset? TIA for any help, so appreciated. This dataset has been a thorn in my side for way too long!! EDITED TO ADD I now think that a tobit model (censored regression) is how I need to model this dataset. It worked well for modeling my fixed effects only. But I still cannot figure out how to make such a model with mixed effects. These datapoints are grouped by site (3 points per site), so it needs to be a random effect. Is there a way to do a censored regression with mixed effects?


title: $\chi^2$ test of significance vs. goodness of fit

text: I have 644 companies out of which 154 went bankrupt. I investigate below if the bankruptcy is related to the sector type:                           Bankrupt     Nothing        total        -----------------------------------     BioTech |     15          110        |  125     Airline |     20          120        |  140     AutoCos |     50          100        |  150     Telecom |     60          40         |  100     Oil&Gas |     9           120        |  129          Hit rate = (15+20+50...) / 125+140+150....      I need to find if the bankruptcy and sectors are inter-dependent. Once that is proven, I have to investigate each sector to check the relationship. I want to check if my following approach is correct: **Test of Independence:** Use (Obs-Exp)^2/Exp on 'bankruptcy' & 'nothing' variables for each group. Sum them up. Expected value for events is HitRate*Total (e.g for BioTech, it is (154/644)x125 for Event. Corresponding expected value for nothing is 125 - (154/644)x125. Use this Chi-Sq and compare it to say 95% probability using **df = 4**. **Goodness of Fit** If test of independence rejects the null hypothesis, I investigate whether each sector has some power. I do exactly what I did above but with **df = 1** for each bin. Q1. Is my second step correct? (i.e. is my goodness of fit test actually a test of independence?) Q2. Is my over all approach statistically correct?


title: Interpreting the sample mean, $\bar y$

text: I have a hypothetical population that contains the values $2,3,4,5,6,7,8,9$. I have to draw a sample of size $4$ from that given population. So I have $\binom{8}{4}=70$ ways to draw the sample. I wrote down all possible samples and computed their sample mean. I got $70$ as the sample mean ($\bar y$). Then I made a frequency table of $\bar y$, showed $\bar y$ is an unbiased estimate of population mean, $\mu_Y$ and drew a histogram of $\bar y$. The histogram shows that the distribution of $\bar y$ is a bell shaped curve. So at the end of this work, my interpretation of the whole task is: > _Though the parent population follows uniform distribution, the sampling > distribution of the sample mean, $\bar y$, tends to Normal distribution._ **Is that all there is to interpreting the whole task? I want to interpret most appropriately so that I can get the whole picture.**


title: Stepwise regression and mediation analysis

text: In my research I have 3 domains of IV, 4 of MV and 8 of DV. To run mediation analysis I made possible chains from the results of stepwise regression and only one chain fulfilled the assumptions of mediation. Then I made the possible chains from the results of correlation and 12 possible chains fulfilled the mediation. Now how it can be possible that a variable is not emerging as a predictor in stepwise regression but it is fulfilling mediation assumption. Or a variable is a significant predictor in stepwise results but it is not significant through enter method of regression?


title: Calculation of power of survival study

text: For an upcoming study of about 200 (rare) cancer cases, we would like to determine the power of detection for a hypothetical marker, present in 10% of cases, i.e. 20 cases with a hazard ratio of 3.0. Cases will be followed up for at least 5 years. We will have the option to validate any identified markers in an independent additional 200 member cohort. In the past, I have used the program PS to run power calculations, which is fairly idiot-proof for a biologist like me. However, the cancer we are studying tends to have an overall survival rate of 60%, so the calculation of a median survival time doesn't really make sense (at least to me). What I wanted to know is whether there are approaches available for the calculation of power which can overcome this limitation, and if so, how can I easily implement them in, ideally, R. I have googled around to find an answer, but haven't met with success so far.


title: Transforming proportion data: when arcsin square root is not enough

text: Is there a (stronger?) alternative to the arcsin square root transformation for percentage/proportion data? In the data set I'm working on at the moment, marked heteroscedasticity remains after I apply this transformation, i.e. the plot of residuals vs. fitted values is still very much rhomboid. Thanks in advance for any suggestions. Edited to respond to comments: the data are investment decisions by experimental participants who may invest 0-100% of an endowment in multiples of 10%. I have also looked at these data using ordinal logistic regression, but would like to see what a valid glm would produce. Plus I could see the answer being useful for future work, as arcsin square root seems to be used as a one-size-fits all solution in my field and I hadn't come across any alternatives being employed.


title: Software for learning Bayesian belief networks

text: I am looking for an easy to use stand alone software that is able to construct Bayesian belief networks out of data. The software should (of course ;-) be free. Can anybody recommend something? Thank you!


title: Reference on how to choose models

text: I wonder if there is any good and easy to understand book on reviewing on when should a particular statistical model (binomial, logistic regression, etc) be used? My background is in statistics, but most of my work is developing MLE estimation algorithm using matrix calculus. I am now starting to focus more on data analysis, which comes in many different ways, and I don't have a good understanding which models are most proper. For example, I know for binary observations I should use logistic regression, but I don't know under which situations I should use a non informative prior (e.g. Jeffreys' prior). Also, since there are so many choice of models, I don't know in a broad way which one to pick. Thanks.


title: Many groups with a single large element and a lot of small ones - How to assess?

text: I have a sample of 100 groups, each with 50 elements. All elements in each group are characterized with a single "size" value. The sizes of all the elements in a group (with 50 elements) are normalized by dividing them with the size of the maximal element (in their group, regardless of the other 49). Hence I have 100 groups, each with at least one element at size 1 (which constitute the maximal element or elements) and the rest are smaller. Observing the data (manually), I notice that in the majority of groups (96 of them) there is a characteristic distribution of sizes : there always exist either a single or two large element (larger than 0.85) and the rest of the sizes are smaller than 0.3. Now, I wonder (and would appreciate suggestions :-) of what would be the right statistical tool to generalize this phenomenon and make claims on the probability that this is the prevailing distribution of elements within any group in the population. Thanks!


title: What to do in logistic regression if you have a huge amount of variables?

text: I am dealing with logistic regression, trying to identify variables which have a causal relationship with a binary response. The way I usually do it is to try variables one by one and visualize the probability of positive outcome curve, and check if it is flat or has a good curve. If it is the latter, then it means there's a causal relationship. I wonder if there is a better routine? Especially if I have a huge number of variables to check, while some of them are not the ones that have a relationship with the observations. What would be the cons if I throw too many variables into logistic regression?


title: Does a median-unbiased estimator minimize mean absolute deviance?

text: This is a follow-up but also a different question of my previous one. I read on Wikipedia that "A median-unbiased estimator minimizes the risk with respect to the absolute-deviation loss function, as observed by Laplace." However, my Monte Carlo simulation results does not support this argument. I assume a sample from a log-normal population, $X_1,X_2,...,X_N \sim \mbox{LN}(\mu,\sigma^2)$,where, $\mu$ and $\sigma$ are the log-mean and log- sd,$\beta = \exp(\mu)=50$ The geometric-mean estimator is a median-unbiased estimator for the population median $\exp(\mu)$, $\hat{\beta}_{\mbox{GM}}= \exp(\hat{\mu})= \exp{(\sum\frac{\log(X_i)}{N})} \sim \mbox{LN}(\mu,\sigma^2/N)$ where, $\mu$ and $\sigma$ are the log-mean and log-sd, $\hat\mu$ and $\hat\sigma$ are the MLEs for $\mu$ and $\sigma$. While a corrected geometric-mean estimator is a mean-unbiased estimator for the population median. $\hat{\beta}_{\mbox{CG}}= \exp(\hat{\mu}-\hat\sigma^2/2N)$ I generate samples of size 5 repeatedly from the LN$(\log(50),\sqrt{\log(1+2^2)})$. The replication number is 10,000. The average absolute deviations I got are 25.14 for the geometric-mean estimator and 22.92 for corrected geometric-mean. Why? BTW, the estimated median absolute deviations are 18.18 for geometric-mean and 18.58 for corrected geometric-mean estimator. The R script I used are here:               #```{r stackexchange}     #' Calculate the geomean to estimate the lognormal median.     #'     #' This function Calculate the geomean to estimate the lognormal     #' median.     #'     #' @param x a vector.     require(plyr)     GM <- function(x){         exp(mean(log(x)))     }     #' Calculate the bias corrected geomean to estimate the lognormal     #' median.     #'     #' This function Calculate the bias corrected geomean using the     #' variance of the log of the samples, i.e., $\hat\sigma^2=1/(n-1)         # \Sigma_i(\Log(X_i)-\hat\mu)^2$     #'     #' @param x a vector.     BCGM <- function(x){     y <- log(x)     exp(mean(y)-var(y)/(2*length(y)))     }     #' Calculate the bias corrected geomean to estimate the lognormal     #' median.     #'     #' This function Calculate the bias corrected geomean using     #' $\hat\sigma^2=1/(n)\Sigma_i(\Log(X_i)-\hat\mu)^2$     #'     #' @param x a vector.     CG <- function(x){     y <- log(x)     exp(mean(y)-var(y)/(2*length(y))*(length(y)-1)/length(y))     }          ############################          simln <- function(n,mu,sigma,CI=FALSE)     {         X <- rlnorm(n,mu,sigma)         Y <- 1/X         gm <- GM(X)         cg <- CG(X)         ##gmk <- log(2)/GM(log(2)*Y) #the same as GM(X)         ##cgk <- log(2)/CG(log(2)*Y)         cgk <- 1/CG(Y)         sm <- median(X)         if(CI==TRUE) ci <- calCI(X)         ##bcgm <- BCGM(X)         ##return(c(gm,cg,bcgm))         if(CI==FALSE) return(c(GM=gm,CG=cg,CGK=cgk,SM=sm)) else return(c(GM=gm,CG=cg,CGK=cgk,CI=ci[3],SM=sm))     }     cv <-2     mcN <-10000     res <- sapply(1:mcN,function(i){simln(n=5,mu=log(50),sigma=sqrt(log(1+cv^2)), CI=FALSE)})     sumres.mad <- apply(res,1,function(x) mean(abs(x-50)))     sumres.medad <- apply(res,1,function(x) median(abs(x-50)))     sumres.mse <- apply(res,1,function(x) mean((x-50)^2))     #```          #```{r eval=FALSE}     #> sumres.mad           GM       CG      CGK       SM      #25.14202 22.91564 29.65724 31.49275      #> sumres.mse           GM       CG      CGK       SM      #1368.209 1031.478 2051.540 2407.218      #```


title: How to choose negative examples to build a training set?

text: I want to use supervised machine learning method to do binary classfication, so a training set is needed. I already have the positive examples. My question is that all my negative examples are randomly selected from a very large dataset which may contain only a few positive examples. Is this approach right to construct the training set? Is any approach to evaluate the training set, especially negative examples?


title: Aggregating pooled regression outputs in different years

text: I performed multiple pooled cross-sectional regressions with the same time- intervals (5years) in different years. I'm wondering now on how to aggregate the different regression outputs. Does it make sense to average the estimated coefficients or do you rather use one simple pooled regression over all $t$s? And assuming taking averages makes sense, how can you get a valid $t$-statistics? Here is an extract of the code I used; thank you for your help.               > head(Grunfeld)       firm year   inv  value capital     1    1 1935 317.6 3078.5     2.8     2    1 1936 391.8 4661.7    52.6     3    1 1937 410.6 5387.1   156.9     4    1 1938 257.7 2792.2   209.2     5    1 1939 330.8 4313.2   203.4     6    1 1940 461.2 4643.9   207.2          library(plm)     data("Grunfeld", package="plm")          # Store each subset regression in myregression     myregression <- list()          # Regression on six-year subsets of Grunfeld     for(t in 1940:1950) {            myregression[[t-1939]] <- lm(inv ~ value + capital,                                    subset(Grunfeld, year<=t & year>=t-5))            # Rename list elements by year range of subset       names(myregression)[[t-1939]] = paste0("Years",t)     }


title: Random or fixed effects in a GLM, or perhaps GEE?

text: I am involved in the analysis of the effect of an intervention at five different hospitals. The intervention aimed to increase the use of a particular equipement when treeating a patient, so the outcome variable is binary (used the equipement/did not use the equipement). The study is not controlled, just pre- and post intervention measurements, with different patients before and after the intervention, so just two independent groups. However, the hospitals are very different with different types of patients. To account for this clustering, I initally considered a logistic regression model with _"equipement use"_ as an dependent variable and _time_ (pre vs post) and _hospital_ as independent fixed effect models, possibly with an interaction term. However, recently I have recieved a suggestion that a GEE would be a prefered choice instead of a logistic regression. I guess this person sees _hospital_ as a random factor, and therefore would like me to consider an population-averaged effect (the population in this case being the hospitals). Does anyone have any thoughts on this? I don't see the advantage of this approach, but I get a feeling that I am missing something crucial. Are there any other suggestions for analysis? The aim of the study is to examine the effect of the intervention, the hospitals per se are not of interest. They are just five convinience sampled hospitals, certainly not a random sample.


title: Minimum cluster size requirements? Minimum sample size requirements for clustered standard errors?

text: I have a sample with little over 100 observations and 50 clusters, one quarter of which have only one observation. Is it correct to calculate clustered standard errors in a linear regression that uses this data or are there some requirements for minimum cluster size and sample size that I am violating? Any references in support of answers provided would be great!


title: What is the statistical justification of interpolation?

text: Suppose that we have two points (the following figure: black circles) and we want to find a value for a third point between them (cross). Indeed we are going to estimate it based on our experimental results, the black points. The simplest case is to draw a line and then find the value (i.e., linear interpolation). If we had supporting points e.g., as brown points in both sides we prefer to get benefit from them and fit a non-linear curve (green curve). The question is that what is the statistical reasoning to mark the red cross as the solution? Why other crosses (e.g., yellow ones) are not answers where they could be? What kind of inference or (?) pushes us to accept the red one? I will develop my original question based on the answers got for this very simple question. ![enter image description here](http://i.stack.imgur.com/5yNZI.png)


title: A/B Testing of a probability, interpreting results of online conversion test

text: So I have a formula to calculate sample size for a Bernoulli trial, given the inputs of the Delta I wish to detect, and the current probability of a success (click-through rate or CTR). $N=21*\theta^2/\delta^2$   95% confidence level and 90% power There are various other formulas one can find that basically mimic this, which are sample size calculators. So I understand everything I have written. I also understand that if set $\delta = 5$%, and that if the A > B by 5% I am only able to conclude that there is a difference (once I meet sample size), but I _can not_ conclude there is a 5% difference. How else can I interpret the results? As in, say I find a 10% difference, how confident can I be that the difference is 10%? Is there anything else I can say besides, "if $\delta \gt 5$% then there is a difference"?


title: What are ways to deal with circular covariates (e.g. with a GAM)?

text: I'm building a model in which several of my covariates live on a "circle", in the sense that they take values in the interval [0,1), and 0=1. I'm wondering about techniques for dealing with this situation. One idea is to represent a circular variable theta as a pair of variables ( sin(theta), cos(theta) ). Any thoughts on this approach or better approaches? I'm specifically using the mgcv package GAMs. Is there a way to tell the model that certain additive pieces should have the same values at the endpoints? Another package? Thanks!


title: How to quickly select important variables from a very large dataset?

text: I have a dataset with about 2,000 binary variables/200,000 rows and I'm trying to predict a single binary dependent variable. My chief goal at this stage isn't getting accuracy of prediction, but rather to identify which of these variables are important predictors. I'd like to drop the number of variables in my final model to about 100. Is there a relatively quick way of getting the most important variables? randomForest seems to be taking a long time. I don't have to use all 200,000 observations, so sampling is an option on the table.


title: Intuition Behind GLM

text: This is a question regarding using logistic regression, and relating it to gaussian distribution or a binomial distribution.               model<-glm(target~ x1, data=data, type='response', family='binomial')     model<-glm(target~ x1, data=data, type='response')  #defaults to gaussian      My understanding of binomial is that it is               theta=chance of success     z=trails ending in success     k=trials ending in failure     (theta^z)*(1-theta)^k      And something Gaussian is               theta = standard deviation     x = success     u = mean     Y = [ 1/σ * sqrt(2π) ] * e -(x - μ)2/2σ2       So I understand how to do GLM with R, I kind of understand what binomial and gaussian means, but I have no understanding of how you relate binomial or gaussian to logistic regression, and how binomial and gaussian are different in this context. Question 1- Can someone explain the intuition behind how "family='binomial'" is used when building a model with GLM? Question 2- Given that the shapes of a binomial distribution and a gaussian distribution look very much the same (they both peak in the middle and gradually go down towards the ends), how does choosing either binomial or guassian lead to different models built from GLM? thanks!!


title: Make Score Value to evaluate social popularity (from fb likes count)

text: I have a list of websites and Facebook likes count for each of them. Count varies from 0 to millions. For each website I want to make up some Score Value from 0 to 10, which would represent sites' social popularity. Any thoughts how to deal with such a problem would be much appreciated. **Update:**   Data summary:   Min. : 1   1st Qu.: 2   Median : 8   Mean : 908   3rd Qu.: 28   Max. :10841643


title: When would time series not work?

text: In my understanding, a time series model is a way to predict the observation. For example, having a time series data points, one possible way to model is using AR model which is $y_{t}$ by $y_{t}= a_{1}x_{t}+a_{2}x_{t-2}+a_{3}x_{t-3}+...+\epsilon_{t} $. There are many applications in real life using time series series model. For example, it seems when people are talking about prediction in finance, it is usually refers to time seres model. However, I wonder why would the time series model would work? It seems really strong assumption on betting the future observation would be a combination of the disturbance, weighted by the coefficients $a_i$. I know that there would be assumptions like stationarity and homoscedasticity (although heteroscedasticity would still keep it unbiased) to ensure the unbiased estimation for the coefficients. However, I guess my problem is at bigger level of asking under what situations the time series model can capture the essence of the target it tries to model.


title: Information theoretic approach to model selection: model diagnostics

text: I've read a few papers that have used the information theoretic to model selection. One of these papers made 36 candidate models. In the methods section of that paper, the authors mentioned nothing about checking the diagnostics (normality, distribution of residuals, outliers etc) of these models. I'm wondering if the researchers bothered to go through the diagnostics of all 36 models. I presume that checking model diagnostics is still something required prior to model selection?


title: CFA Fit Measure Thresholds for New Scales

text: A well established questionnaire measuring _interpersonal_ relations is modified to measure _interorganizational_ relations. 11 parameters of the relation are determined by 4 items each. Items are removed until the best possible $\alpha$ per scale is achieved (.58 - .81). This results in some scales with 3 items, some scale with 4 items. A confirmatory factor analysis is applied to calculate the estimates of the items _and_ verify the internal structure of the questionnaire. As expected with some $\alpha$ as low as .58 and .61 the CFA does not fit very well. **Since it is "kind of" a new scale, is it ok to use the predictions for further modelling and report them?** * * * The fit indices are (bad):                 Number of observations                           251            Estimator                                         ML      Robust       Minimum Function Test Statistic             2475.330    2029.585       Degrees of freedom                               854         854       P-value (Chi-square)                           0.000       0.000       Scaling correction factor                                  1.220         for the Satorra-Bentler correction          Model test baseline model:            Minimum Function Test Statistic             5718.372    4583.881       Degrees of freedom                               946         946       P-value                                        0.000       0.000          Full model versus baseline model:            Comparative Fit Index (CFI)                    0.660       0.677       Tucker-Lewis Index (TLI)                       0.624       0.642          Loglikelihood and Information Criteria:            Loglikelihood user model (H0)             -18173.075  -18173.075       Loglikelihood unrestricted model (H1)     -16935.410  -16935.410            Number of free parameters                        180         180       Akaike (AIC)                               36706.150   36706.150       Bayesian (BIC)                             37340.731   37340.731       Sample-size adjusted Bayesian (BIC)        36770.110   36770.110          Root Mean Square Error of Approximation:            RMSEA                                          0.087       0.074       90 Percent Confidence Interval          0.083  0.091       0.070  0.078       P-value RMSEA <= 0.05                          0.000       0.000          Standardized Root Mean Square Residual:            SRMR                                           0.146       0.146      * * * ## Edit: The real question is probably not only whether to continue with estimates based of a bad fit, but rather whether the scale is new or not. * * * ## Edit2: Barrett (2007) recommends to not rely on "approximate fit indices" at all. > The criterion used for "fit" is actually an abstract concept in the majority > of SEM models. It is clearly not predictive accuracy. In fact whether models > "approximately fit" with an RMSEA of 0.05 or 0.07 is a literally meaningless > scientific statement. * * * # Bibliography Barrett, P. Structural equation modelling: Adjudging model fit. _Personality and Individual Differences_ , **2007** , _42_ , 815-824 [PDF]


title: Why at all consider sampling without replacement in a practical application?

text: Sampling with replacement has two advantages over sampling without replacement as I see it: 1) You don't need to worry about the finite population correction. 2) There is a chance that elements from the population are drawn multiple times - then you can recycle the measurements and save time. Of course from an academic POV one has to investigate both methods. But from a practical POV I don't see why one would consider sampling without replacement given the advantages of with replacement. But I am a beginner in statistics so there might be plenty of good reasons why without replacement might be the superior choice - at least for specific use cases. Please, unconfuse me!


title: Can adjusted-$R^2$ compare models across different samples?

text: I was reading this journal article about the decrease in the relevance of the income statement: Basically what the author does is this: (1) For each year, he takes that year's (say, 2010) income statement information from a set of companies as the $X$ values (independent variable), and the 1-year return (2011 share price/2010 share price) of the companies' stocks as the $Y$ values (dependent variable). He generates a regression model and calculates the adjusted $R^2$ for this model. (2) He does the same in part (1) for 30 years (3) He plots a graph of the adjusted $R^2$ for each year over time, and generates a regression model for it. It turns out that adjusted $R^2$ has decreased by about 10% to 7%. Thus, he concludes that the income statement has lost relevance over time. I've seen adjusted $R^2$ used to compare different models using the same sample data. But I've never seen it used in such a manner before because the data used to generate each model is different. Any advice would be much appreciated.


title: How to explain the logistic regression model objects from the rms::lrm() function in R?

text: Recently I was trying to do logistic regression using the `rms::lrm()` function. But I had some trouble understanding the model objects from the function. Here is the example from the package:               #dataset     n            <- 1000    # define sample size     set.seed(17)            # so can reproduce the results     treat        <- factor(sample(c('a','b','c'), n,TRUE))     num.diseases <- sample(0:4, n,TRUE)     age          <- rnorm(n, 50, 10)     cholesterol  <- rnorm(n, 200, 25)     weight       <- rnorm(n, 150, 20)     sex          <- factor(sample(c('female','male'), n,TRUE))     L            <- .1*(num.diseases-2) + .045*(age-50) +                     (log(cholesterol - 10)-5.2)*(-2*(treat=='a') +                     3.5*(treat=='b')+2*(treat=='c'))     y            <- ifelse(runif(n) < plogis(L), 1, 0)     #fit model     g            <- lrm(y ~ treat*rcs(age))          > g          Logistic Regression Model          lrm(formula = y ~ treat * rcs(age))                                Model Likelihood     Discrimination    Rank Discrim.                                  Ratio Test            Indexes          Indexes            Obs          1000    LR chi2      76.77    R2       0.099    C       0.656          0            478    d.f.            14    g        0.665    Dxy     0.312          1            522    Pr(> chi2) <0.0001    gr       1.945    gamma   0.314         max |deriv| 3e-06                          gp       0.156    tau-a   0.156                                                    Brier    0.231         > anova(g)                     Wald Statistics          Response: y            Factor                                     Chi-Square d.f. P           treat  (Factor+Higher Order Factors)        5.62      10   0.8462       All Interactions                           1.30       8   0.9956      age  (Factor+Higher Order Factors)         65.99      12   <.0001       All Interactions                           1.30       8   0.9956       Nonlinear (Factor+Higher Order Factors)    2.23       9   0.9872      treat * age  (Factor+Higher Order Factors)  1.30       8   0.9956       Nonlinear                                  0.99       6   0.9858       Nonlinear Interaction : f(A,B) vs. AB      0.99       6   0.9858      TOTAL NONLINEAR                             2.23       9   0.9872      TOTAL NONLINEAR + INTERACTION               2.57      11   0.9953      TOTAL                                      69.06      14   <.0001      **Here are my questions:**   For the object `g`,   * What does the `max |deriv| 3e-06` mean?    * What do the Discrimination and Rand Discrim. Indexes suggest?  For the `anova(g)` object,   * What's the `Factor +Higher Order Factors` for the treat?    * Why there are two `all interactions`? How to explain the nonlinear parts?


title: Comparing distributions of generalization performance

text: Say that I have two learning methods for a **classification** problem, $A$ and $B$, and that I estimate their generalization performance with something like repeated cross validation or bootstrapping. From this process I get a **distribution of scores** $P_A$ and $P_B$ for each method across these repetitions (e.g. the distribution of ROC AUC values for each model). Looking at these distributions, it could be that $\mu_A \ge \mu_B$ but that $\sigma_A \ge \sigma_B$ (i.e. the expected generalization performance of $A$ is higher than $B$, but that there is more uncertainty about this estimation). I think this is called the **bias-variance dilemma** in regression. What **mathematical methods** can I use to compare $P_A$ and $P_B$ and eventually make an informed decision about which model to use? **Note:** For the sake of simplicity, I am referring to two methods $A$ and $B$ here, but I am interested in methods that can be used to compare the distribution of scores of ~1000 learning methods (e.g. from a grid search) and eventually make a final decision about which model to use.


title: What's the difference between a confidence interval and a credible interval?

text: Joris and Srikant's exchange here got me wondering (again) if my internal explanations for the the difference between confidence intervals and credible intervals were the correct ones. How you would explain the difference?


title: Whether or not to use standardization in this example

text: here is my scenario: People record data using sensors, the data is uploaded to a server. Then the users can view their data on a web page. Here is what is in my dataset: - I know how often each user uploaded a record to the server - I know if and how often that user viewed his records on the website I want to know how often the users viewed their records on the website. Now since not all users have the same life-style and available time to make records and upload them: Would it make sense to standadize the data? If you think yes how would I do that? Wikipedia tells me the standard scores are calculated by z = (x - mean)/standard deviation. So x is my number of record views on the web page, mean is the mean of the views of all users and standard deviation is ... well the standard deviation of all users' views of their data? Any help is appreciated. If you think I am completely on the wrong path please tell me =)


title: Whether to use an offset in a Poisson regression when predicting total career goals scored by hockey players

text: I've got a question concerning wheter or not to use an offset. Assume a very easy model, where you want to describe the (overall)number of goals in hockey. So you have goals, number of games played and a dummy variable "striker" which is equal to 1 if the player is a striker and 0 otherwise. So which of the following models is correctly specified?   1. goals=games+striker , or   2. goals=offset(games)+striker Again, the goals are overall goals and the number of games are overall games for a single player. For example there could be a player picked up who has 50 goals in 100 games and another player who has 20 goals in 50 games and so on. What am I supposed to do when I'd like to estimate the number of goals? Is it really necessary to use an offset here? References:   * See this previous question discussing when to use offsets in Poisson regression in general.


title: Calculating the information contained in a message

text: Lets say i want to calculate the information content of a particular message.What apart from the message itself has to be taken into account in doing so, and what data would i need to collect to perform my action?


title: Rationale of using AUC?

text: Especially in the computer-science oriented side of the machine learning literature, AUC (area under the receiver operator characteristic curve) is a popular criterion for evaluating classifiers. What are the justifications for using the AUC? E.g. is there a particular loss function for which the optimal decision is the classifier with the best AUC?


title: Goodness of fit in GLMs

text: I am searching for a good criterion to measure the "goodness of fit" in generalized linear models. To make clear: I am not searching for a criterion which gives me an answer to the question "does overdispersion occur?". What do you think about Nagelkerke's pseudo R-squared? Any thought would be appreciated!


title: Memory requirements of $k$-means clustering

text: Can anyone tell me the factors that affect the memory requirements of $k$-means clustering with a bit of explanation?


title: Sum of normally distributed variables is a normally distributed variable?

text: Consider two Wiener processes: $$ \begin{aligned} X_a &\sim\mathcal N(0,a) \\ X_{a-b} &\sim\mathcal N(0,a-b) \end{aligned} $$ How do I show that: $$ X_a - X_{a-b} \sim\mathcal N(.,.) $$ That is, how do I show that $X_a - X_{a-b}$ is Gaussian? * * * I know that $E(X_a - X_{a-b}) = 0$ and $$ \begin{aligned} \text{Var}(X_a - X_{a-b}) &= E(X_a^2) - E(X_{a-b}^2) \\ &= \text{min}(a,a) -\text{min}(a-b,a-b) \\ &= a - a - b \\ &= b \end{aligned} $$ * * * Note that both $X_a$ and $X_{a-b}$ have the following properties:   * They are normally distributed with mean $0$.   * Their sample paths are continuous.   * Its covariance function is $\rho(s,t) := \text{Cov}(X_s,X_t)=E(X_sX_t)=\text{min}(s,t)$.


title: Verification of a regression model

text: I need some guidance related to regression model verification using validation data. I am new to R-tool & statistics and trying my best to learn. I did search on internet too but I couldn't get a final answer to my questions. Actually I have a lot of questions, I may try my best to explain the problems: I am experimenting with network packets and R-tool. I have captured some packets from a network using a custom made packet sniffer in java. The sniffer will capture some packets and save the information of packet header like: tcp window size, tcp sequence numbers, date-time, ip header length, ip time to live etc... in a csv file. Also the sniffer will add category number to each csv file so that we can know which packet belongs to which category. I created 9 different categories saved in 9 different csv files. Now I extracted 1000 observation from each of the csv files and created a data set named "alldata". Then I created training data set and validation data set from "alldata" data set. Now I want to perform linear regression, logistic regression, decision tree analysis, cluster analysis etc on this "alldata" data set. So my plan was to use training data set to create models and then later use validation data set to verify my models. Category will be my target variable in any case. I want to predict the category from other independent variables.   1. My first confusion is that after I created scatter plot of category with other independent variables and I don't see any linear relationship between them. Moreover I even don't know what relation exists between category and independent variables. From scatter plots it seems to me that there is no specific relation between category and other independent variables(except date_time it is bit linear to category). Am I doing the correct interpretation ? Here are some of the plots: plot 1 plot 2   2. I think doing linear regression won't make any sense now after having a look at scatter plots. Is this correct assumption?   3. Although I tried to do make some regression models with training data set, but the R-square values for all the models is quite low (for example like 0.00019, 0.0035, 0.018 etc. ) So can I assume that these models are not good due to very low r-square vales?   4. As logistic regression is used when we have target variables having only two values 1 or 0, or some probabilities between 0.0 - 1.0. This means performing logistic regression is not possible for this type of data set. Is my assumption true?   5. My main question was how to verify a model created with training data set by using validation data set? Please let me know the commands and the procedure. Please let me know if I am doing this in wrong way or if you can suggest me a better way to do this whole work. I think if someone could please clear my doubts then I may ask further more questions. If you don't understand my problem we can discuss in more detail I look forward for your replies. Thank you! * * * @Wayne Hello thanks for the reply, but the thing is for each category I have almost same range of values of independent variables like(tcpheader, ipttl, iplen). For example iptype is only having two values 6 and 17. So most of the categories are having iptype value of 6 & 17. So it is also same is for tcpheader, tcp sequence number, tcp acknowledgement number etc. I don't think there is any way to distinguish a particular packet based on these independent variables. Only the independent variable that can be helpful is time. But when I created a model with time it had good r-squared value but the regression line equation doesn't predict category with any value of date_time. I don't understand this behaviour. Thanks.


title: Using z-Score for multiple comparisons

text: I have data collected from three groups and I would like to calculate the z-score using the grand mean and variance of the combined three groups. Would it be possible to then compare the z-scores of each subgroup and thereby check if they all 3 groups are taken from one population? Would this method allow me to circumvent multiple comparisons correction?


title: Any actual applications based on "dirichlet distribution servers as a conjugate prior for multinomial distribution"

text: Conjugate prior is useful and beautiful in theory, for instance, Dirichlet distribution can serves as the conjugate prior for multinomial distribution. But I don't find any actual applications based on that. Can you provide some?


title: Adaboost for neural networks. Is it still worth it?

text: I have a question about Adaboost and neural networks. Given the recent development in neural networks (dropout, maxout, or rectified linear units) is there a significant benefit of performing Adaboost in order to increase accuracy? I am asking the question also having in mind that in order for Adaboost to perform well you need to train multiple neural networks which is not very time efficient and if the improvement in accuracy is not substantial then maybe it is not worth it. Thanks!


title: How to chose optimal number of lags and inputs?

text: I'm using Genetic Algorithms to do inputs selection in a time series problem. The issue is that the number of possible inputs is very large (100 possible inputs + inputs' lags) and I don't know a procedure to determine the minimum number of lags to consider. My procedure was to run a correlation analysis ex ante and leave de top 70 inputs in terms of correlation with the desired output. The problem with this procedure is that, as I increase the number of possible lags in consideration, more "ridiculous" inputs (i.e lag 140) are selected...


title: Generalized least squares with insignificant predictor variable

text: Suppose I have fitted an standard linear regression mode $Y=\beta_0+\beta_1X_1+\beta_2X_2+\beta_3X_3+\epsilon$. Based on the ACF plot or PACF plot of the residuals for this regreesion model, I found that the $\epsilon$'s are serially auto-correlated that can be modeled as an AR(1) i.e. $\epsilon=\phi\epsilon_{t-1}+\nu_t$, where $\nu_t$ are white noise. I decided to fit the above model with Generalized least squares approach with the following code               library(nlme)     model.gls=gls(Y~X_1+X_2+X_3,data=dat, correlation=corARMA(p=1), method="ML")      A part of output in R is as follows:               Correlation Structure: AR(1)      Formula: ~1       Parameter estimate(s):           Phi      0.8067566           Coefficients:                     Value  Std.Error    t-value p-value     (Intercept) -3.172158 0.01504180 -210.88950  0.0000     x1          -3.890222 0.14141947  -27.50839  0.0000     x2           0.089368 0.11667430    0.76596  0.4456     x3          -0.155348 0.10190694   -1.52441  0.1308      Clearly $X_2$ and $X_3$ are not significant based on the t-test. Q1: Is it a correct approach to drop $X_2$ and $X_3$ and refit a new Generalized least squares model with the **same** correlation structure (i.e. AR(1))? In this case the code should be like:               model.gls=gls(Y~X_1,data=dat, correlation=corARMA(p=1), method="ML")      I am asking this because the AR(1) has been selected while all $X_1$, $X_2$ and $X_3$ have been included in the linear regression model. So dropping two terms may change the correlation structure (i.e AR(1)) to another ARMA model. If the approach in Q1 is not correct, how can we proceed and handle the non- significant terms, $X_2$ and $X_3$ in the fitted Generalized least squares? Thank you.


title: How to elegantly group subjects?

text: I want to conduct an experiment to find out whether the CO2-production of soil changes with high amounts of rainfall. Therefore, I want to group 8 soil samples into two different groups. One group is control, one group is treatment. Due to time, space and money restrictions, I can only have four replicates for each group. I have some pre-experiment CO2-production rates of these 8 samples that I would like to use to group the subjects. The two groups should be as similar as possible before the treatment starts (equal mean and variance). How can I programmatically determine the "best" division into groups? Data of pre-test CO2 production are below. Any direct help or linkage to other websites/posts is greatly appreciated. Please let me know if and how I can improve the question. Some google-fu was to no avail.               #### Data: y1-y8 are the soil samples, the values are CO2 procution rates     y1 <- 10     y2 <- 20     y3 <- 22     y4 <- 30     y5 <- 15     y6 <- 12     y7 <- 28     y8 <- 26          group1.1 <- c(y1,y2,y3,y4)      group1.2 <- c(y5,y6,y7,y8)          group2.1 <- c(y4,y7,y8,y3)     group2.2 <- c(y1,y2,y5,y6)          boxplot(group1.1, group1.2,las=1) # not too bad, but I am sure it can be done better...but how?     boxplot(group2.1, group2.2,las=1) # bad...


title: Problem with Bioconductor/limma on two-color factorial design: Coefficients not estimable

text: I'm running into a problem where my coefficients are not estimable when using Bioconductor/limma on a two-color factorial design for a microarray analysis. I have microarray data I downloaded from Array Express using the ArrayExpress function in the ArrayExpress package. I managed to convert the NChannelSet object into an RGList object that I can use in limma using a tip I found in the archives (https://stat.ethz.ch/pipermail/bioconductor/2009-September/029705.html). After some bg subtraction and normalization, I have an MAList object with 32 arrays (red and green) by 34,944 probes. The experiment consists of two cell types (WT or KO), by 4 treatments (control, CD70, CD80, CD70+CD80), at four different time points (2, 4, 8, 14 hours), with a dye-swap. I'm mostly interested in genes that are differentially expressed when WT cells are hit with CD70+CD80 versus CD80 alone. I'm trying to create a targets file so I can use modelMatrix to create my design matrix. But I have a multifactorial design using a two color dye-swap design, and I'm not sure how to specify this in the targets file. The limma manual has information about factorial designs, but no examples for two-color experiments. I have the following factors: 1: Celltype: WT or KO 2: Treatment: control (x), CD70, CD80, or CD70+CD80 (CD7080). 3: Timepoints: 2, 4, 8, and 14 hours. I tried collapsing all of factors and levels into a single string in my targets file:               > targets            array          Cy3            Cy5     1  array3329    WT.CD70.2         WT.x.2     2  array2675  KO.CD7080.8    WT.CD7080.8     3  array2242  WT.CD7080.2    KO.CD7080.2     4  array3328       WT.x.2      WT.CD70.2     5  array3310       WT.x.8      WT.CD70.8     6  array2246  KO.CD7080.2    WT.CD7080.2     7  array3337  WT.CD7080.4      WT.CD80.4     8  array3323 WT.CD7080.14     WT.CD80.14     9  array2673    KO.CD70.8      WT.CD70.8     10 array1938 WT.CD7080.14   KO.CD7080.14     11 array2240    WT.CD70.2      KO.CD70.2     12 array3336    WT.CD80.4    WT.CD7080.4     13 array3322   WT.CD80.14 WT.CD7080.14.2     14 array2674  WT.CD7080.8    KO.CD7080.8     15 array3321   WT.CD70.14        WT.x.14     16 array2241    KO.CD70.2      WT.CD70.2     17 array2597  KO.CD7080.4    WT.CD7080.4     18 array3313  WT.CD7080.8      WT.CD80.8     19 array1939 KO.CD7080.14 WT.CD7080.14.2     20 array3335    WT.CD70.4         WT.x.4     21 array2672    WT.CD70.8      KO.CD70.8     22 array3320      WT.x.14     WT.CD70.14     23 array3334       WT.x.4      WT.CD70.4     24 array3311    WT.CD70.8         WT.x.8     25 array3331  WT.CD7080.2      WT.CD80.2     26 array1941   KO.CD70.14   WT.CD70.14.2     27 array2588    WT.CD70.4      KO.CD70.4     28 array2596  WT.CD7080.4    KO.CD7080.4     29 array3330    WT.CD80.2    WT.CD7080.2     30 array3312    WT.CD80.8    WT.CD7080.8     31 array1940 WT.CD70.14.2     KO.CD70.14     32 array2593    KO.CD70.4      WT.CD70.4      I then created a design matrix, where the reference group is WT cells, untreated, 2-hour timepoint.               design <- modelMatrix(targets, ref="WT.x.2")      This creates a 32x25 design matrix, but when I fit the model, I get a "Coefficients not estimable" error/warning:               > fit <- lmFit(d, design)     Coefficients not estimable: WT.CD70.14.2 WT.CD80.14 WT.CD80.2 WT.CD80.4 WT.CD80.8 WT.x.14 WT.x.4 WT.x.8      Warning message:     Partial NA coefficients for 34944 probe(s)       Any help with how I can specify this multifactorial time-course design in a two-channel dye-swap experiment would be greatly appreciated. As I said, I'm most interested in the conditions on arrays 8-9, 13-14, 19, 26, and 30-31, where I'm looking at WT cells treated with CD70 and CD80 versus CD80 alone. Thanks in advance for any help!


title: Choosing between chi-squared / logistic regression vs difference of mean tests for studying bankruptices

text: You'd expect more companies to go bankrupt if they have low cash balances. Of course, bankruptcies can also happen even when you have lot of cash on hand (to terminate labor contracts, etc.). The distributions are shown below: ![enter image description here](http://i.stack.imgur.com/zn4Iw.jpg) A `chi-squared` test shows an association of bankruptcies with cash-on-hand. I then check the results using `logistic regression` which actually gives me negative coefficients for low cash balances! (Using LR is fine as I have ~2000 positives for ~150k observations). I'm thinking about a simple `difference of means test` to show that cash balances for bankruptcies differ than those for healthy companies. Do you suggest any other statistical test that I can use instead?


title: What does T statistics of information coefficient indicate?

text: Hi I am looking for a clear explanation of T statistics concept. Especially in quantitative equity portfolio management context, what does T statistics of information coefficient indicate? How can I tell if it is significant, and if it is, any useful information I could derive from that?


title: Differences in coefficients

text: Suppose I want to see whether $z$ is a confounder for a model with $y$ the outcome variable and $x$ the predictor. If I adjust for $z$, and the adjusted coefficient of $x$ changes versus the unadjusted coefficient of $x$, does it matter by how much it changes? If the difference between the unadjusted and adjusted coefficients for $x$ is very small, should I still take this into account?


title: Calculating inter rater reliability where raters and ratees only partially overlap

text: I'm new to inter rater reliability calculations. I have 5 developers and 6 raters. The raters ranked the first 3 of these developers based on some criteria like the following example: For example,                        dev1  dev2  dev3  dev4   dev5     ranker1   1     2     3     ranker2   3     1     2     ranker3         3     1     2     ranker4               1     3      2     ranker5   2     1           2     ranker6         3     1            2      The ratings measure the "level of expertise": "1" the best expertise, "2" less expertise, and "3" lowest expertise My question is what the best algorithm to calculate the inter rater reliability among these raters?


title: PCA to decorrelate variables

text: I have 2 variables that I want to decorrelate. I was told I can use PCA to do so. I did PCA on the data and got all the parameters. Now how do I get the new set of transformed data that no longer correlates with the second variable? I want to use this set for further analysis. thanks,


title: Efficient numerical methods to estimate ARMA models?

text: I mean exact likelihood based estimation instead of these LS methods. There are more general nonlinear optimzation methods, but in terms of performance, are there any specific methods for this type of problems? thanks.


title: A best measure for speaker recognition

text: I have a set $E_{1}$, with a finite cardinality $n$ of rectangular matrices which contains the useful MFCC coefficients generated from $n$ speech signals. Similary I have a set $E_{2}$ of same cardinality as that of $E_{1}$ which is a collection of vectors of finite dimension containing the LPC of the same set of speech signal which was used to form $E_{1}$. Now $ D=\{ E_{1},E_{2} \}$ forms the database for the speaker recognition system. When a test signal is given, MFCC $M_{i}$ and LPC $L_{i}$ are generated and the closest members $M_{j} \in E_{1}$ for MFCC and $L_{j} \in E_{2}$ for LPC are found using a distance function $d$. It is not necessary that $M_{j}$ and $L_{j}$ are the exact members of $E_{1}$ and $E_{2}$ respectively. It depends on the acoustic environment during the test phase.   * What is the distance function used in literature?    * If its $L_{2}$ norm, is there any other better measure which is "sensitive" so I can reduce the possibility of misclassification?


title: R only alternatives to BUGS

text: I am following a course on Bayesian statistics using BUGS and R. Now, I already know BUGS, it's great but I am not really fond of using a separate program rather than just R. I have read that there are a lot of new Bayesian packages in R. Is there a list or reference on which packages there are for Bayesian statistics and what these do? And, is there an R package alternative for the flexibility of BUGS?


title: Is there an R package with a pretty function that can deal effectively with outliers?

text: One of the data sets I deal with is quite strange. The datawarehouse I downloaded the data from has a lot 999999999 values in one of the variables. Apparently the computer system on which the datawarehouse sits on does not support storing of null values. So they use 999999999 as the "null" value. Now if I just run `pretty` in R on the variable, it gives non-sensical ranges.   * Is there a package with a version of `pretty` that can deal with outliers by putting them in the range of say (100,High) ?


title: Statistical power and minimum sample size for ANOVA with likert scale as dependent variable

text: I'm trying to compute the minimum sample size for a psychometric test based on 7 point Likert scales. I'd like to run ANOVA on each scale to look for differences between groups. Most online survey sample size calculators seem to be designed for polls, e.g. Yes/No, Agree/Disagree. They take as input population size, a confidence interval and a proportion (50% Yes/50% no) and then return the required sample size. Most statistical books suggest using power tests (such as R's power.t.test), which take as input a minimum effect size, alpha, beta and a statistical test and then return the required sample size. For my purposes power tests seems to make the most sense, but what has me concerned is that none of them take into account the population size, which seems like it ought to have at least some effect on the outcome. So my question is, what is the correct calculation to use in my specific survey situation and more generally what is the connection between power tests and these online survey sample size calculators, does population size matter in some way, perhaps helping to capture the notion of representative sample?


title: Simulation of maximum likelihood ratio test to test two poisson random variables

text: I have two random poisson variables $x_1$ and $x_2$ with value 10 and 25 respectively. I am interested to use likelihood ratio test to test the null hypothesis: $\lambda_1=\lambda_2$, versus alernate hypthesis $\lambda_1$ not equal to $\lambda_2$. I want to use simulation to calculate power and alpha values. I would want to do it in R so any reference to R codes will be appreciated. Thanks in advance


title: Three-way repeated measures ANOVA - what to do next?

text: I ran a 4(A)x4(B)x2(C) ANOVA in SPSS. I got significant main effects of A, B and C, and significant interaction effects of AxB and BxC. I got no significant interaction effects of AxC and AxBxC. My questions are: What do I do now? What post-hocs do I carry out? Should I make plots for AxB and BxC disregarding the third factor?


title: How to find +/- uncertainty with a least squares regression

text: I have a set of data points with an uncertainty on each point. From these data points I can fit a line, who's which the slope is a significant value. How do I use the information I have to get an estimate on the $ \sigma $ of my final value? All I can find online is information on $ R^2$ and $ \chi ^2 $. Bonus points if there's an easy way to do it in python/numpy/scipy.


title: Statistics test for concentration of extract effecting germination

text: I am undertaking biology course work and am a little stuck. My experimental hypothesis is: There is a significant correlation between increasing tomato extract concentrations and reduced seed germination. Which statistical test would be best?


title: Modeling nearest event through multivariate random variable with condition

text: There are two time sequences in my system, one of them represents IN events, and the other OUT events. Each IN event would be released by the nearest following OUT event or reaching the deadline. I want to get the pdf of the interval between the time the IN event is triggered until the time it is released (OUT event occurs or censored at deadline). As the following figure shows, $t_a$, $t_r$, $t_r^'$ represent the IN event time, the nearest OUT event time and last OUT event time prior to $t_r$. The nearest OUT event can be represented as $t_r>t_a$ and $t_r^' < t_a$. I extend the situation to statistics, so that $t_a$ is a sample of random variable $T_a \sim U(0,T)$. $t_r$ obey $T_r \sim U(0,2T)$. The interval of two OUT event, $t_r-t_r^'$, should be treated as $T_i \sim Exp(\lambda)$. ![enter image description here](http://i.stack.imgur.com/RfuzT.png) How to model this complex multivariate random variable and get its pdf? Furthermore, deadline of event $t_a$ is $2T$, which means that if there is no $t_r$ between $(t_a, 2T)$, it would be released at time $2T$. This condition could introduce more difficulty of modeling.


title: What are different measures that say something about a given dataset?

text: Given a data set (having _n_ instances and _m_ features), what are the different measures that give some insight about some properties of the data set? In other words, if there are two such data sets, are there any measures that talk about similarity/difference between the data sets? Thank you.


title: Whether to include individual items or the composite variable in structural equational modelling?

text: I made performing structural equation modelling. Some of my latent variables include a lot of items (e.g., 20 or more items). Alternatively, I could just include in the model a composite variable based on the items that I have computed myself. **When performing SEM, is it better to include all the items or just the composite variable?**


title: Are there linear multiple regression methods resilient to missing data?

text: Are there any multiple linear regression methods or packages that are resilient to occasional missing values? I have no prior view on imputing the missing values based on the nature of the data, and I would like to avoid discarding rows that have NAs. Although I am not performing a panel regression, my data is arranged as panel data: Date, Identifier for individual in population, Characteristic 1, Characteristic 2 , ... , Objective function value.


title: What does the index variable k define in the Lasso regularization function

text: In the Lasso L1 regularization, from where comes the value of the variable $k$ in the second part of the function? Why isn't it $n$, too? $$L(\beta) = \sum_{i=1}^n (y_i - \phi(x_i)^T \cdot \beta)^2 + \lambda \sum_{j = 1}^k l(\beta_j)$$ $\beta$ is the parameter vector, $y$ the output vector and $x$ the input vector.


title: Missing data and discriminant analysis

text: I have a question. A company conducted a study in two countries using the same battery of questions. However, one question was left out in country B by mistake (the questions ranged in a scale from 1 to 5 and they were sort of: how important is for you to have a product that is reliable, meets my needs, etc.). After the battery of questions was complete, discriminant analysis were used to compute membership to an attitudinal group (initial groups came from a previous study in a third country, and from this study, the coefficients of the discriminant were obtained for each question). When people realised the mistake, they proposed three solutions:   * to take the mode of the other variables (there were 20 in total) and substitute the value of the missing question with that one for each person asked (this assumes lack of independence)   * to make the coefficients of the discriminant analysis to zero for that particular missing question   * to generate random numbers for this particular question for each person. I wonder if any of these solutions actually make sense?


title: Standardization of compositional data in PCA versus using real data

text: I have a question about conducting a PCA between variables that are measured in different units. I understand the importance of using a correlation matrix versus a covariance matrix to minimize variance. The data I'm working with is not normally distributed and has not been transformed in other tests. For example, there are three variables A, B, and C, and 20 observations, where 10 observations are measured using 1 set of units, and the other 10 observations are measured using another set of units*. The values between the units are quite different in in value and variance (expected). The data is not normal in either units and has not been transformed. The measurements using the first set of units is 2 to 3 orders of magnitude higher than those measured using the other units (expected). I have conducted a PCA using a correlation matrix and interpreted results. However a non- statistician recommended I `standardize' the measured data, such that I'm using ratio or fractions for all the observations for each of the variables: _Variable A/(sum of all 3 variables) and so on and so forth for Variables B and C._ However, a PCA using the contributing fraction of each variable is different from PCA using the measured value in different units leading to different eigenvalues and eigenvectors, thus leading to two different scientific interpretations. Beyond using the appropriate association matrix, Is this "standardization" step valid and or necessary from a statistical perspective? **Update:** Should PCA be done on compositional data?


title: Comparing two dependent intra-class correlations (ICCs)

text: We are measuring two waves from a biological signal. As an example consider the ECG $P$ wave and $T$ wave. From a single trace, the amplitude of both of these waves is measured. We want to compare the inter-session reliability of the two waves; doing a second measurement a week later. So in this example we have $P_1$, $P_2$, $T_1$ and $T_2$. We expect, and have found, correlations between ($P_1$, $P_2$) ($P_1$, $T_1$) and all of the combinations. $ICC_{2,1}$ for the $P$ wave ($P_1$, $P_2$) was $0.8$ and for the $T$ wave was $0.6$. Confidence intervals for each ICC can be computed by bootstrap but how can we compare the test-retest reliabilities of the two waves?


title: What are the sampling distributions of higher moments of the normal distribution?

text: Let $X_i$ be independent, normally distributed random variables, for $1\leq i\leq N$. What is the distribution of $Y_m=\frac1 N \sum_{i=1}^N X_i^m$? Every high school student knows part of the answer. The mean of $Y^m$ is the $m$'th moment of the normal distribution, and the variance of $Y^1$ is $1/N$. I'm interested in the width of the distribution of $Y^m$. How does its variance scale with the sample size $N$? (Is variance a useful measure?) I'm especially interested in the large sample limit, with m a small integer. Knowing where to look it up would be helpful.


title: p value correlation for random vectors in R

text: I want to calculate correlation coefficient between random vectors and use this for hypothesis testing. The code               bs <- function(qtS, rvNum = 1000)     {       eb <- qtS * rvNum       randVect <- matrix(rnorm(eb), nrow=qtS, ncol = rvNum)       randVect <- randVect/sqrt(sum(randVect * randVect))       randVectCor <- cor(randVect)       randVectCor <- randVectCor[col(randVectCor) < row(randVectCor)]       plot(density(randVectCor), main = "resampled data density distribution", type = "l", lwd  = 2, xlim = c(min(randVectCor), max(randVectCor)))       list(min = min(randVectCor), max = max(randVectCor), randVectCor = randVectCor)     }      works well but with one peculiar thing that I argue a lot with my colleague. Namely, if the qtS is small (i.e. 5) then the max and min are higer than 0.99 which in turn leads to rejection of the vast majority of possible observed correlation values. But with increasing qtS min and max values get lower and with i.e. 50 this value is 0.6. Is this OK, I mean the fact that possible correlations between vectors with fewer elements are orderly larger than between vectors with greater number of elements?


title: Hidden Markov model alternative that allows dependence on previous n states

text: I am working on a prototype framework. Basically I need to generate a model or profile for each individual's lifestyle based on some sensor data about him/her, such as GPS, motions, heart rate, surrounding environment readings, temperature etc. The proposed model or profile is a knowledge representation of an individual's lifestyle pattern. Maybe a graph with probabilities. I am thinking to use Hidden Markov Model to implement this. As the states in HMM can be Working, Sleeping, Leisure, Sport and etc. Observations can be a set of various sensor data. My understanding of HMM is that next state S(t) is only depends on previous one state S(t-1). However in reality, a person's activity might depends on previous n states. Is it still a good idea to use HMM? Or should I use some other more appropriate models? I have seen some work on second order and multiple order of Markov Chains, does it also apply HMM? I really appreciate if you can give me a detailed explanation. Thanks!!


title: Do the results show an omitted variable bias? How sizeable is the bias?

text: The question is more involved on how to calculate the omitted variable bias. We were given data and are supposed to use SAS to run regression models. First, how do you know if results suggests there was an omitted variable bias? Second, how do you estimate the different components of Bias($\tilde{\beta_{1}}$)=$E(\hat{\beta_{1}})-\beta_{1}=\beta_{2} \tilde{\delta_{1}}$ since you don't know the true values of $\beta_{1}$ and $\beta_{2}$?


title: Details of Bayesian linear regression

text: In Bayesian Linear Regression, we have a data set with {$x_{i}$, $t_{i}$} where $x_{i}$ are input vectors and $t_{i}$ are their resulting observations. We want to find a vector $\bf{w}$ in order to maximize the posterior distribution $p({\bf{w}}|{\bf{t}})$ and hopefully make good prediction for future $x_{i}$'s. In this Bayesian formalism, we multiply a prior $p({\bf{w}}|\alpha)$ by the likelihood function $p({\bf{t}}|{\bf{w}})$ to obtain $p({\bf{w}}|{\bf{t}})$ and then use this distribution as a prior and repeat the process again. However, I have a question regarding the details of this procedure. For example, suppose we don't have any data points at first and rely exclusively on the prior, and as long as we can infer the hyperparameter $\alpha$, we are going to be able to do that. In the next step, we receive a data point $\{x_{1}$, $t_{1}\}$, and we calculate the likelihood $p({\bf{t_{1}}}|{\bf{w}})$ which involves a Gaussian with mean ${\bf{w}}^{T}\phi({x_{1}})$ and variance $\beta^{-1}$. Using a model such as $w_{0}+w_{1}x$, this amounts to this: $$p({\bf{t}}|{\bf{w}}) \propto \exp{\{\beta(t_{1} - (w_{0}+w_{1}x_{1}))^{2}\}}$$ so if that is correct, we know everything except $w_{0}$ and $w_{1}$. Do I get a probability for each combination of $w_{0}$ and $w_{1}$? In any case, I think this is a function of $\bf{w}$ that we can multiply by the prior to obtain a posterior. Did I understand anything wrong about this procedure of $\text{posterior}\propto\text{likelihood}\times\text{prior}$? If so, I would appreciate any corrections. **UPDATE** : I added a picture from Bishop's book "Pattern Recognition and Machine Learning" in which this process is described. The first column is the likelihood, second column is the prior/posterior and third column describes the input vectors in blue and samples of $\bf{w}$ taken from the posterior in the equation $w_{0}+w_{1}x$. ![enter image description here](http://i.stack.imgur.com/5Zayv.png) Thanks


title: How to test difference between two items rated on a Likert scale by same participants?

text: I have conducted a survey with multiple 5 point Likert scale questions ranging from 1 = very unimportant to 5 = very important. I used frequency distributions, and stacked bar charts to show which percentage of respondents rated the items as either important or very important. There are some great differences between the importance placed on the different items: e.g. one item was only rated (very) important by 6 % of the respondents while another one was rated (very) important by 70 % of the respondents. To me, this difference seems very significant. **Is there a statistical test that I can run to determine whether the difference between those two t-point likert scales items is statistically significant?**


title: How to write the derivative of Expected Shortfall

text: I am trying to implement a gradient of Expected Shortfall/CVaR in $R$ According to Kalkbrenner, the derivative ($\Pi_{\alpha}^{ES}$) of ES is defined as: $$ \Pi_{\alpha}^{ES} = E_{\mathbb{Q}_X}(X_i) = \dfrac{( \int X_i \ast \mathbb{1}_{\{X > q_{\alpha}(X)\}} d\mathbb{P} + \beta_X \int X_i \ast \mathbb{1}_{\{X = q_{\alpha}(X)\}} d\mathbb{P})}{1- \alpha} $$ I'm basically failing to grasp the formula here. Could anyone explain how I would approach this, I have a feeling it's not that hard in practice. * * * **edit:** I did some crunching on my own.. since $X$ is the main portfolio stochastic variable, $X_i$ is the i-th portfolio. And if the last part $\beta_X \int...$ can be assumed to be $=0$ if no value exist exactly _on_ the $\alpha$-quantile (I dare to assume that) What is left is more or less the simple Marginal Expected Shortfall, which I was trying to replace in the first case. I will probably stick with MES unless someone here proves me wrong


title: Using conjoint analysis for price elasticity

text: I have been asked to use conjoint analysis to find the price elasticity of a certain category of products. I have never used conjoint analysis before and have been reading about it for a while. I have to confess I'm not delighted with the methodology and I wanted to know if anyone has any experience with this. Can anyone endorse this approach or mention some pros - cons resulting from their experience? Even more so, what would people actually use to do something like this? Not sure if this is relevant but I would be working in R.


title: What are the implications of a perfect fit model?

text: I perform logistic regression with a relatively small dataset (N=65), using 12 parameters (11 variables, one constant, no interactions), which results in a perfectly fitting model (in SPSS). I have a few questions about this:   1. Is it correct to say that the model saturated, or is that only the case if there are less data points than parameters?   2. Is there a difference between a saturated model and an overparameterized model?   3. What should I do, if anything, with the fact that this does model result in a perfect fit? I perform this test as a robustness check, iterations of the model that leave out some of the parameters corroborate the main findings. I don't necessarily _need_ coefficient estimates of this model, because I'm mostly interested in the signs and I obtained those from a previous test (same regression, but differently defined variables). Thank you in advance.


title: Help with understanding this covariance setup

text: I have been reading a paper that formulates the problem of image registration as a generative model and I have been having a lot of trouble understanding some concepts and I was wondering if someone would be kind enough to shed some light on a particular matter that has me stumped for a few days now. The problem aims to find a non-linear geometric transformation that aligns the two images geometrically (i.e. find a dense correspondence between the image points). This is an ill-posed problem and needs some regularisation. So, if the geometric transformation be parameterised by 'w', then we need define some sort of a prior on the distribution of these transformation parameters. In this work, it is done by specifying the prior as a multivariate normal distribution with 0 mean and a covariance matrix as follows: $$ P(w) = N(w; 0, (\lambda K)^{-1}) $$ where $K$ is the inverse covariance matrix scaled by $\lambda$. Now, one way to promote plausible transformations in image registration is penalising some measure of the "energy" of the transformation. A commonly used one is something called the bending energy which is a sum of the directional second derivatives and can be written down in a matrix form and the bending energy can be calculated as: $$ E = w^T M w $$ where M would be the bending energy matrix. Now in this work, what the authors have done is taken this bending energy matrix (which has a sparse form) and used that as the inverse covariance matrix i.e. $$ P(w) = N(w; 0, (\lambda M)^{-1}) $$ Now, I am not sure how this can be used in place of the cover. matrix. From what I see as the definition of the covariance matrix, an entry in the cell (i, j) is given by $E[(w_i-\mu_i)(w_j-\mu_j)]$. I guess my broad question is how does one translate any such sort of commonly used regularisation schemes like penalising second derivatives of transformations into a valid covariance or inverse covariance matrix. How should one go about specifying this matrix? Sorry for the long, convoluted question. I would appreciate any help you can give me. Thanks, Luca


title: Derivation for the area formula for a sample data taken from a population that is normally distributed?

text: I need to understand where this formula comes from and what starting assumptions are required to get there: $ f_i = \frac{i-0.375}{n+0.25}$ where $i$ is the index (the position of the data value in the ordered list) and $n$ is the number of observations. The expected proportion of observations less than or equal to the $i$th data is $f_i$. $f_i$ represents the area to the left of $z_i$, our z-score. Context: A normal probability plot is a graph that plots observed data vs. "normal scores." A "normal score" is the expected z-score of the data value assuming that the random variable $X$ is normal.


title: Tukey honestly significant difference

text: Can I use Tukey's HSD not as a post hoc test? just instead of a serious of T-TESTS? I am checking 8 variables. Also 5 of my variables as dichotomous and 3 have 3 levels Thanks


title: Pattern Clustering

text: I have data of daily time spent in studying for around 2000 students. I need to make clusters based on the pattern of hour spent in a week, not by average hour spent.For example some students are studying for 5 hours a day during weekdays & in weekend it is 7 hours,for some students 2 hours in weekdays & weekend it is 4 hours & for some students 0 hour in weekend & 6-8 hours in weekdays.I want to cluster these patterns along with some other macro economic variables. Can you please suggest me the logic for pattern clustering or link to read about pattern clustering? Thanks in advance


title: Comparing different step functions

text: I measured the appearance of certain events over time in different conditions and I would like to know whether there is a difference between the happening of these events. I thought of plotting the events as cumulative curves, such as:               # The times at which events occurr     times.1 <- c(0, 5, 15, 20, 25, 33, 40, 50)     times.2 <- c(0, 10, 20, 21, 35, 40, 56)     # The number of events     num.1 <- c(0, 10, 30, 22, 11, 3, 8, 10)     num.2 <- c(0, 5, 3, 20, 10, 8, 5)     plot(times.1, cumsum(num.1), t="s", col="red", lwd=2,          xlab="Time", ylab="Number of events", las=1)     points(times.2, cumsum(num.2), t="s", col="blue", lwd=2)      This gives me this nice graph: ![cumulative graph](http://i.stack.imgur.com/zvnXe.png) How would analyse whether there are differences between the two conditions, knowing that I have 5-10 measurements for each group? Also, I would like to represent an average curve for each condition, but the events do not always happen at the same time, so I find myself a little bit stuck as I cannot think of some sensible way of "averaging" discrete non uniform data. I am interested both in changes in the number of events and in the timing (e.g.: a certain condition may reach the same number of events but slower than the other)


title: Interpretation of MSE (mean square error) and ME (mean error)

text: I'd like to produce forecasts by considering different scenarios. I'm using Mean Error (ME), where the error $=$ forecast $-$ demand, and Mean Square Error (MSE) to evaluate the results. For the scenarios that bias (ME) is negative the MSE is very high, how can I interpret these results? I know that the MSE$=$variance of forecast error + bias$^2$, so for these scenarios, we have a low bias but MSE is high, so it means the variance of forecast error is high. ![enter image description here](http://i.stack.imgur.com/Kpbd2.jpg) ![enter image description here](http://i.stack.imgur.com/tnxBi.jpg)


title: Correct approach to statistical analysis of voting results

text: Consider data consisting of voting results. Votes can be either individual votes of jury members or aggregated votes by certain communities (professional of geographical). In the end there is a set of data where each row has three fields: Voter, Candidate (the one who collects votes) and Place (or rank, voters put candidates on different places). There are several competitions, so candidates are measured multiple times. Now we want ot analyze if some voters are biased, e.g. if they keep on giving some candidates certain places no matter if those candidates deserve it. I am not really sure what's the proper way of achieving this. Here is what I am trying to figure out:   1. In case votes are collected from communities, voting results are based on raw data of different sizes: one community can have 100K members while the other one has only 1000. Can we still assume that "places" represent ordinal data, not categorical?   2. What is the recommended technique to evaluate if a voter is biased? I can think of a simple measurement: an average difference between a candidate's place and the place he is given by a certain voter. The more the difference, the more the candidate is overrated by this voter. But I am not sure this is right thing to do because with such approach a single very low mark can balance multiple overrated high marks.   3. I've read about chi-square and Fisher's tests but not sure if it's applicable here. AFAIK they are applicable to categorical data. Any help is appreciated.


title: Robust logistic regression vs logistic regression

text: It's been a while since I've thought about or used a robust logistic regression model. However, I ran a few logits yesterday and realized that my probability curve was being affected by some 'extreme' values, and particularly low ones. However, when I went to run a robust logit model, I got the same results as I did in my logit model. Under what circumstances should a robust logit produce different results from a traditional logit model? (in terms of coefficients) R Code:               > library(Design)     > ddist<- datadist(dlmydat)     > options(datadist='ddist')     > me = lrm(factor(dlstatus) ~ dlour_bid, data=dlmydat)     > me          Logistic Regression Model          lrm(formula = factor(dlstatus) ~ dlour_bid, data = dlmydat)               Frequencies of Responses       1   2      906 154                  Obs  Max Deriv Model L.R.       d.f.          P          C        Dxy      Gamma      Tau-a         R2      Brier            1060      3e-05     170.11          1          0       0.81      0.619      0.621      0.154      0.263      0.105                     Coef      S.E.      Wald Z P     Intercept -5.233549 0.3731235 -14.03 0     dlour_bid  0.005367 0.0004925  10.90 0          > library(car)     > dlmod = glm(factor(dlstatus) ~ dlour_bid, data=dlmydat, family=binomial(link="logit"))     > summary(dlmod)          Call:     glm(formula = factor(dlstatus) ~ dlour_bid, family = binomial(link = "logit"),          data = dlmydat)          Deviance Residuals:          Min       1Q   Median       3Q      Max       -1.2345  -0.5687  -0.3059  -0.1739   2.6999            Coefficients:                   Estimate Std. Error z value Pr(>|z|)         (Intercept) -5.2335492  0.3731235  -14.03   <2e-16 ***     dlour_bid    0.0053667  0.0004925   10.90   <2e-16 ***     ---     Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1           (Dispersion parameter for binomial family taken to be 1)              Null deviance: 878.61  on 1059  degrees of freedom     Residual deviance: 708.50  on 1058  degrees of freedom     AIC: 712.5          Number of Fisher Scoring iterations: 6


title: Glasserian Grounded Theory Building Blocks

text: On reading the FAQ I was in hope this question pertaining a particular aspect of grounded theory would fall under data analysis. I appreciate if you would point me to which other community I should refer this to otherwise. My question is of a novice on the field of Grounded Theory. I have seen and read few papers that distinguish what is and what is not, and the difference between Glausserian and Straussian Grounded Theory. In an attempt to perform Grounded Theory from a Glausserian (or Classical as it is called) point of view, I've been trying to find what are the meaning of their building blocks. In their website the following building blocks are defined: `indicator, incident, category, code, concept` My question is thus of asking their meaning. Although one could believe that checking a book from them to be the simplest way, such literature, at least from the point of view from experts in my domain (Software Engineering), has been even considered to be mystical.


title: Use pc scores to check mulitvariate normality

text: How can you check multivariate normality, using the scores from the PCA? Or what can we expect about the scores if the data is multivariate normal distributed?


title: RFM & customer lifetime value modeling in R

text: Can anybody tell me how to do recency, frequency & monetary value (RFM) modeling & customer value modeling in R? Also, can somebody refer me some literature on it?


title: Power analysis for matched poisson variables

text: I have two matched groups consisting of 22 observations each for a Poisson distributed variable. I am interested in evaluating the power to detect a 20% decrease in the outcome for the second group. I have no idea of the correlation between matched groups and hence will have to assume independence. I am planning on evaluating the power using a simulation study and have the following ideas in mind. 1) Define rate0 and rate1 = rate0-(0.20*rate0) 2) Simulate X1~Poisson(22,rate0) and X2~Poisson(22,rate) 3) Compute pairwise differences X1-X2 and conduct a paired Wilcoxon test or a paired ttest. 4) Conduct simulation nsim times and compute the times pvalue < 0.05 I would appreciate some feedback/suggestions regarding step 3. Would anyone have any other ideas?


title: Getting sense of and testing multiple correlations in two-block structure

text: I am exploring a multidimensional data set. The variables can be mainly divided into two groups: A and B. I have calculated the correlation (Spearman's rho) matrix of this data set, and found that although the correlation between a variable from group A and one from group B is generally low, there are nonetheless some interesting patterns. For instance, variable A1 in group A has a (small yet) _positive_ correlation with each variable in B, A2 has a (small yet) _negative_ correlation with each variable in B, etc. Are these patterns potentially useful? Is there a way to formally measure them?


title: How to calculate advertising copy wearout using kalman filter?

text: I have to calculate half life of an advertisement using Kalman filter in R. The paper 'Estimating the Half-life of Advertisements' (Naik, 1999[1]) provides the base but am unable to understand how exactly the estimation is being done and since a lot of parameters that go in as the input to the kalman filters are unknown i am unable to move forward with my project. Any suggestion as in how to go forward? [1]: Naik, P.A. (1999), "Estimating the Half-life of Advertisements," Marketing Letters 10:3, 351-362


title: Subspace clustering with random transformation

text: One approach for clustering a high dimensional dataset is to use linear transformation, and the most common approaches are PCA and random projection (where random projection arises from the Johnson-Lindenstrauss Lemma). I was wondering why we can't use other random transformation s like when our transformation matrix $R$ was drawn from a uniform distribution?


title: How to fit a cumulative distribution over two proportions at different values

text: I am trying to figure how how I would fit a distribution to the data described below. Imagine I run an experiment where in condition 1 subjects choose between option A costing \$50 and option B costing \$100. In condition 2 subjects chose between option A costing \$50 dollars and option B costing \$80. So imagine that I get 20% preferring option B in condition 1, and 40% preferring option B in condition 2. To me that says that 20% of the population values B over A by at least \$50 and 40% of the population values option B over A by at least \$30. Let's say I assume a normal distribution around how much people value B over A. It seems like if I know that 20% value the difference over \$50 and 40% value it over \$30 I should be able estimate the total distribution of this difference in the population. I was also thinking that maybe I need one more condition to properly estimate. Is this correct? If so, what would be the right tool to estimate that distribution? Thank you!


title: Computing probability 1 player win given both follow normal distribution with same standard deviation

text: I am given that score ("logscore") of 2 players follow a normal distribution. They both have the same standard deviation of 0.205 but the mean varies. Player A has mean logscore that exceeds B by 0.15. I need to compute the probability A wins B. How do I approach this problem? ![enter image description here](http://i.stack.imgur.com/hJME8.png) * * * I think the its also very similar to the below question (e). How do I do such comparisons on 2 normal random variables? ![enter image description here](http://i.stack.imgur.com/20zrE.png)