title: Auto-thumbnail generator model based on interval measures for tube sites

text: I think video streaming services e.g. tube like sites could become the next playground for data scientists, particularly I was quite interested in play rate enhacing and want to do some research about this topic. It came to me some ideas on how user-systen information could be analyzed in the way of building something like an auto-thumbnail generator model that selects the highest density frame of the clip, this would involve working with datasets where each measure looks more like an interval. Let's say a clip starts on second 0 and user1 decides to skip directly to second 120 after watching only the first 15 seconds and then almost inmediatly jumps backs to the beginning - let's say second 10. The data would look something like this: 0:15 120:121 10:X I wonder if there exists any kind of framework or technique to work with this sort of intervalar measures, hopefully my explanation makes sense.


title: How would you explain logit coefficients to someone with no statistical background

text: Similar question here. The question I highlighted above provides an overview of how generalized linear models work. However, I find people often want more and ask about how the coefficients were calculated, and how the coefficients are then turned into odds-ratios or predicted probabilities. So if I were explaining logit models to someone with no experience with statistical modeling, how would you explain the process of going from data to coefficients, and then on to odds ratios and/or predicted probabilities with as limited math as possible?   1. How does the model take data and produce specific coefficients?   2. How do we get useful information from these coefficients? It's difficult to avoid math on this one, as it's important to point out how the probability curve or odds ratios are calculated. My question deals more with 1. than 2.


title: Multinomial regression with categorical choice and predictos and factor analysis

text: I have a non-ordinal categorical dependent variable with 3 choice outcomes and 20 ordinal categorical predictors and want to do a multinomial logistic regression. However, I want to reduce the predictors to just a few variables with factor analysis. SPSS gives regression factor scores as variables for each of my cases. Do I replace my original 20 predictors with these 5 regression factor score variables, or do the factors need to be on another format before using them as predictors?


title: Dirichlet process vs. Pitman-Yor and power-law - how do they relate?

text: The Dirichlet process has, on expectation, $\alpha \log n$ "restaurants", when borrowing the Chinese restaurant process terminology, with $\alpha$ being the concentration parameter. Pitman-Yor, on the other hand, has $\alpha n^d$ "tables", with $\alpha$ being the concentration parameter, and $d$ being the discount parameter. How does that relate to power-law? I understand PY has power-law behavior, but DP does not. From Wikipedia: "The parameters governing the Pitman–Yor process are: 0 ≤ d < 1 a discount parameter, a strength parameter θ > −d and a base distribution G0 over a probability space X. When d = 0, it becomes the Dirichlet process. The discount parameter gives the Pitman–Yor process more flexibility over tail behavior than the Dirichlet process, which has exponential tails. This makes Pitman–Yor process useful for modeling data with power-law tails " This, unfortunately, is not clear enough for me.


title: PCA on out-of-sample data

text: We know that the projection matrix learned by PCA can be applied to out-of- sample data points to get their low-dimensional embedding. However, how reliable are these embeddings expected to be, as compared to the embedding obtained from PCA with these out-of-sample points combined with the original data? Consider this hypothetical pseudo out-of-sample setting: Let's say I have 1000 data points and I want to do PCA on them. Can I instead do PCA on maybe just 500 of them (in order to have some computational savings) and then use the learned projection matrix to get the embedding of the rest of the points as well (by treating them as out-of-sample data)?


title: What does it mean in terms of regression if residuals are not white noise?

text: I need help in answering this one, it is an exam question.


title: Interpreting the neural network output in R?

text: I managed to create neural network of my data. But I am not so sure about the interpretation of the R output. I used following command to create neural network:               > net=nnet(formula = category~iplen+date_time, size=0,skip=T,lineout=T)     # weights:  3     initial  value 136242.000000      final  value 136242.000000      converged      Then I used following command to see the output:                   > summary(net)     a 2-0-1 network with 3 weights     options were - skip-layer connections       b->o i1->o i2->o       0.64 -0.46  0.15      So from the above output Can I can conclude the following diagram of neural network?: ![net](http://i.stack.imgur.com/sg4Ji.jpg) Second question is how can I know how useful this diagram is? I mean I wanted to find the category number(target variable) from the independent variables. so now how can I decide if this network really helped me to predict the category(target variable)? What is the final output or how to find that? Can I conclude the following output from the above n-network? :               category= -0.46(iplen)+0.15(date_time)+0.64      Regards,


title: Should I reverse score items before running reliability analyses (item-total correlation) and factor analysis?

text: I want to assess item-total correlations on a 19-item questionnaire (some of the questions are meant to be reverse-scored). My question is:   * Do I reverse score the items PRIOR to calculating the item-total correlations (in order to eliminate any variables that do not correlate with the total at >.40)?    * Additionally, should the items be reverse-scored prior to running a factor analysis?


title: Classification and regression trees (cart)

text: In `matlab`, classregtree can be used to implement classification and regression trees (CART) you can find this in the documentation however it's not clear what methods are used for either classification or regression, 3 methods exist: **Bagging** decision trees, an early ensemble method, builds multiple decision trees by repeatedly resampling training data with replacement, and voting the trees for a consensus prediction. **A Random Forest** classifier uses a number of decision trees, in order to improve the classification rate. **Boosted Trees** can be used for regression-type and classification-type problems. Rotation forest - in which every decision tree is trained by first applying principal component analysis (PCA) on a random subset of the input features. I'm tempted to say random forest but would like to clarify, does anyone know which method matlab implements for classification?


title: How can you select one of two individuals at random using a biased coin

text: Suppose you are given a biased coin for which the probability of getting a head is $p (0<p<1)$.Discuss how you will select one of two individuals at random using the biased coin.


title: degrees of freedom in ARIMA model

text: I estimated parameters of an ARIMA model by arima function in R and I calculeted t ratios or t statistics for each parameter. Now I want to find p values for t test, what is my degree of freedom? is it the number of my observations?


title: Is there a fuzzy distinction between discrete and categorical variables?

text: Are there some variables where it can be difficult to assess whether they are purely categorical? **For example: ** To assess the level of pain that a patient is in, you use a scale between 1 and 10.I wouldn't consider pain score to be a categorical variable here, just discrete; this is because it is a measurement (even if it is subjective). Conversely if I have something such as tumour stage, it is still a discrete variable between say 1 and 4, however I am more tempted to say these are categories. The thing is, increasing value in the Tumour stage does actually have meaning: more likely to have aggressive disease. For categories such as hair colour, or gender, there is no numerical relationship between the categories! The only thing I can really see is that in discrete variables, each step increase has the same meaning (thus increase of 1-2, 5-6 and 9-10) is the same. Whereas with tumour stage va the step sizes don't mean the same thing : 1-2 is different than 2-3 and 3-4 because the tumour aggression is probably not linear. This has implications when running logistic regression models, where you could treat categorical variables differently (names you code them differently).


title: Interpretation of intercept when adjusting for baseline response

text: I am working on a data set with two time points and I am interested in the change from baseline. There are no differing treatment groups. First, I realize that this is most simply made using a paired t-test. However, I would like to adjust for baseline response. My simple question is: Am I right to interpret the intercept term in the following model as the mean baseline adjusted change? y= score at time2 x= score at time1 $\alpha$ = intercept $$y - mean(x)= \alpha + \beta (x- mean(x))$$


title: Is there a general method for simulating data from a formula or analysis available?

text: De novo simulation of data from an experimental design data frame.   With a focus on R (though other language solution would be great). In designing an experiment or a survey, simulating data and conducting an analysis on this simulated data can provide terrific insight into advantages and weaknesses of the design. Such an approach can also be essential to the understanding and proper use of statistical tests. However, this process tends to be somewhat tedious and many are led to skip past this important step in an experiment or survey. Statistical models and test contain most of the information required to simulate the data (including an assumption or an explicit statement of the distribution). Given an an analysis model (and its associated assumptions eg. normality and balance), the levels of a factor and a measure of significance (such as p-value), I would like to obtain simulated data (ideally with a generalized function akin to print(), predict(), simulate()). Is such a generalized simulation framework possible? If so, is such a framework currently available? Example, I would like a function, such as:                sim(aov(response~factor1+factor2*factor3),               p.values=list(factor1=0.05,                             factor2=0.05,                             factor3=0.50,                             factor2:factor3=0.05),               levels=list(factor1=1:10,                           factor2=c("A", "B", "C"),                           factor3=c("A", "B", "C")))      ie, a generalized version of:               sim.lm<-function(){     library(DoE.base)     design<-fac.design(nlevels=c(10,3,3),                        factor.names=c("factor1", "factor2", "factor3"),                        replications=3,                        randomize=F)          response<-with(design, as.numeric(factor1)+                           as.numeric(factor2)+                           as.numeric(factor3)+                           as.numeric(factor2)*as.numeric(factor3)+                           rnorm(length(factor1)))          simulation<-data.frame(design, response)}      OR               sim(glm(response~factor1+factor2*factor3, family=poisson),              p.values=list(factor1=0.05,                            factor2=0.05,                            factor3=0.50,                            factor2:factor3=0.05),              levels=list(factor1=1:10,                          factor2=c("A", "B", "C"),                          factor3=c("A", "B", "C")))      OR                 library(lme4)       sim(lmer(response~factor1+factor2 + (factor2|factor3)),                F_value=list(factor1=50,                             factor2=50),                levels=list(factor1=1:10,                            factor2=c("A", "B", "C"),                            factor3=c("A", "B", "C")))      that would create a complete corresponding data.frame potential examples of specific functions (please edit at will)   \- arima.sim function exist to create a data.frame of the factor levels, without the modelled response:   eg. conf.design   http://cran.r-project.org/web/views/ExperimentalDesign.html


title: Interaction terms

text: Suppose you have variables $X_{1}, \dots, X_{100}$ (can be discrete or continuous). If you want to test all possible interactions, should you include them all in 1 model? We are doing linear regression here.


title: Estimating prices

text: I am currently trying to estimate what the prices for exchanging dollars to Lempiras (Honduran money) and price of fuel will be on June 30, 2012. The data is very limited, Honduras isn't particularly good at keeping records of things, so all I have is what the prices have been over past months and years, how can I use this data to get a valid estimation?


title: Cross validation for transformed and untransformed outcome

text: Suppose that I am interested in selecting from one of the following model spaces: $$\begin{align*} y &= \beta_0 + \beta_1 x_i + \epsilon_i \, \, \text{or} \\\ \log (y_i) &= \gamma_0 + \gamma_1 x_i + \eta_i. \end{align*}$$ I intend to use cross validation for making my selection. If I calculate the prediction errors on the validation set, the first model will give errors in $y$, while the second will give errors in $\log(y_i)$. Obviously, it doesn't make sense to compare these. Instead, for the second model, suppose that I estimate $$\begin{equation*} \text{E}[y_i \mid x_i] = \exp \left\\{\hat{\beta}_0 + \hat{\beta}_1 x_i \right\\} \times \exp \left\\{ \hat{\sigma}^2 /2 \right\\}, \end{equation*}$$ using an assumption that $\epsilon_i \mid x_i \sim N(0, \sigma^2)$ (alternatively, I could just use a Poisson model with robust standard errors/oversidsperion to get the same result; see here for a nice explanation). Then, I compare this predicted $y$ to the observed validation $y$ to get my prediction errors. Is this procedure sensible? Will it generally have good properties?


title: Sum of absolute values of T random variables

text: Where X is a r.v. following a symmetric T distribution with 0 mean and tail parameter $\alpha$. I am looking for the distribution of the n-summed variable $ \sum_{1 \leq i \leq n}|x_i|$. $Y=|X|$ has for PDF $\frac{2 \left(\frac{\alpha }{\alpha +y^2}\right)^{\frac{\alpha +1}{2}}}{\sqrt{\alpha } B\left(\frac{\alpha }{2},\frac{1}{2}\right)}$, $y \geq 0 $. I managed to get the characteristic function $C(t)$ but could not invert the convolution, that is, $C(t)^n$. Thank you for the help.


title: Showing MGF of Poisson converges to MGF of N(0,1)

text: I'm trying to finish up a proof of the CLT for the Poisson distribution, but am having some trouble evaluating a limit. I've shown that the moment generating function for the standardized Poisson is $$M_X (t) = \exp \left \\{-t\sqrt{\lambda}-\lambda+\lambda e^{t/\sqrt{\lambda}} \right \\}$$ (this is the form of the MGF which was asked for us to show) But I'm having some trouble evaluating $$\lim_{\lambda \rightarrow \infty} \exp \left \\{-t\sqrt{\lambda}-\lambda+\lambda e^{t/\sqrt{\lambda}} \right \\} = e^{t^2/2}$$ Any suggestions for my mathematically-challenged mind?


title: How can I test Control variables within a mediation analysis?

text: In my Masters thesis I do a mediation analysis with Multiple Regression Analysis. So I make use of the Baron & Kenny method for testing mediation. I also want to control for some variables. For example I want to know if gender has an effect on the mediation effect. What I did was to split the data file(women/ men) (I'm using SPSS 21) and then run the mediation analysis again separately for the two groups. However, if I want to control for age, splitting the file makes no sense. Is there any analytic procedure to test the effect age or any other non categorical variable has on the mediation effect? And if this should be the case, how can I run this analysis in SPSS? Thank you very much for your help.


title: Calculating the PMF of roll-and-keep dice (XdZkY)

text: I've been trying to use the PMF shown in an older question, found here, and sourced from here, but I've run into a problem. When attempting to calculate N, I'm finding that B (using s+i+j-x-ry = B) rapidly becomes negative when keeping more than 1 die (y>1). This, in turn, results in attempting to distribute a negative number of "balls". This is a bit out of my depth, since I'm not very experienced with combinatorics, but as far as I know I understand all the operations I'm trying to perform. What am I missing in using this? I'm not sure what to make of the modifications from the original Physics Forum version to the version in the question. I've been using the original Physics Forum one, where they differ, but at the moment my problem lies in a part that seems to be identical.


title: Law of Large numbers and central limit theorem

text: I have come across a bit of a dilemma in my exam revision, this is not a topic I am particularly strong with so help is most appreciated: Assume that $X_1,X_2,...$ is an i.i.d sequence, where $E(X^2)=5$ and $E(X)=0$, Using central limit theorem, find the distribution of $$ \frac{\sum_{i=1}^nX_i}{\sqrt{5n}} $$ Here I have let $ S_n= \sum_{i=1}^nX_i $ and we know that $E(X)=0$, then we have: $$ Yn=\frac{\sum_{i=1}^nX_i - nE(x))}{σ\sqrt{n}} = \frac{\sum_{i=1}^nX_i}{\sqrt{5n}} -> N(0,1) $$ Next I am asked to apply the strong law of large numbers to $$ \frac{\sum_{i=1}^nX_i^2}{n} $$ to show that it converges almost surely and compute its finite limit I have had a go at this by letting $ S_n= \sum_{i=1}^nX_i^2 $ and $ E(S_n)= \sum_{i=1}^nE(X_i^2) $ and plugging this into to the strong law of large numbers formula ( $\frac{S_n-E(S_n)}{n}-> 0$ a.s.) to find $ \frac{\sum_{i=1}^nX_i^2}{n} $ -> $ E(X_i^2) $ a.s. However I do not feel as though I am showing that $ \frac{\sum_{i=1}^nX_i^2}{n} $ is converging almost surely using this method? Any ideas would be appreciated Also lastly I have to find, by using central limit theorem, the distribution of $$ \frac{\sum_{i=1}^nXi}{\sqrt{\sum_{i=1}^nX_i^2}} $$ I am lost on this one completely, although I am making the assumption as we are using central limit theorem with it them it will be N(0,1) but as I say that is just a guess, Help/hints/explanations would be really appreciated


title: Research study design inclusive of mediators

text: At workplace, team members often have disagreements with each other. On the one hand disagreements facilitate group discussion, and stimulate creative ideas but on the other hand, disagreements damage team cohesiveness, and decrease team performance. As a researcher, you are interested in this workplace phenomenon, and want to conduct an empirical study on such topic. Please describe your research design by answering the questions below.   1. Clarify what are the conceptual variables (constructs) in your study?   2. Clarify what are the operational variables in your study?   3. Briefly formulate your research hypotheses.   4. How are you going to measure or manipulate the relevant variables?    5. Which statistics methods would you like to use in building a causal relationship?    6. Are there any special points which you would like to address in your research design?


title: multiple hypothesis testing for Poisson distribution

text: I have count data (observed) and expected counts calculated by regression analysis. I want perform multiple hypothesis test to check which elements show maximum differences between observed and expected. I am using R and I ran this comand,               p.values=apply(matrix,1,function(x){t.test(x[1],x[2],alternative="g")$p.value})      It gives me error saying `not enough x observation` -it is obvious that I don't have any replicates. But the paper which I am following (they have done similar kind of study) says that: > To detect hotspots, we performed a one-tailed test (alpha<0.001) using > Poisson distribution and the expected counts I don't understand how did they do it. I can detect the elements using fold-change between observed and expected, but here they say they have done hypothesis test.


title: Statistical significance of long term trend?

text: I have data sets of long term climatology wind speed data on daily basis for 1950-2000. I have calculated the mean for every year and calculated trend pattern, it shows some positive trend of 0.05 m/s/year. I have been asked do the statistical significance test for this trend analysis. How could I achieve it?


title: E-step of the stochastic approximation EM

text: I am reading the paper: Convergence of a stochastic approximation version of the EM algorithm to implement this algorithm for a probability model I already have. In p. 3, the paper summarises the algorithm as follows. ![enter image description here](http://i.stack.imgur.com/lxL66.png) I am stuck at the E (or S) step in this algorithm. In a typical EM setup, one maximizes the integral $$Q(\theta)= \int \log p(x,y |\theta) p(x|y,\theta) dx$$ In here, this integral is estimated via the samples simulated from the posterior. I understand this. But, what I do not understand is: Do authors propose to 'update' the cost function? If so, how can we maximize the new $Q$? May be I am missing a very obvious thing and can not see how to implement this algorithm (I do not understand updating the 'cost'). Thanks in advance!


title: Stratified sample strategy - sample size

text: We have conducted a survey and collected over 30,000 surveys, from school children, ages between 13-18 years old, on a national level across 18 cities, spanning 3 regions (east, south and west), of a country, let's say X. We calculated the sample size as 2396 (based on country population, and percentage of people aged 15-24). The data will be analysed on age and gender. I'd like to know if my following strategy is correct for calculating the stratified surveys, before analyzing the response rate. Apologies, if I'm using incorrect terms, I do not have a background in statics and new to all this.   1. Firstly, I calculated the total percentage of surveys collected from each city.    2. Second, I then calculated above data based on sample size. To get sample for each city.   3. Third, I calculated the original data collected split into gender then into age groups, and proportionally downsized to sample requirements.  Please see attached example below. ![My calculations for one city](http://i.stack.imgur.com/qjQmC.jpg)


title: Why are log transforms used in measuring dose-time association?

text: Explain why log(dose) and log(time) are used to measure the association between dose and time instead of the actual values of dose and time.


title: How to interpret R-squared in multiple regression with more sets of dummies and continuous variables?

text: I have a problem with a multiple regression I performed:   * model without constant term;   * one dependent continuous variable;   * first set of dummies: derived from 2 continuous variables, I used the median value of them as a threshold to derive two binary variables; from these two binaries, I derived 4 dummies, one for each combination (10, 01, 00, 11);   * second set of dummies: 3 dummies derived from one categorical variable;   * two continuous variables. This model has a r-squared value of 98% (and similar adjusted r squared): I think it is too high, but I don't know how to interpret it correctly and assess its eventual validity; I know that r squared tend to increase with the number of explanatory variables, but I don't know if the number of dummies has an influence in its value and validity as an indicator of a good regression. Moreover, this model present high VIF values, indicating collinearity: are these measures still valid or not? I have to say I have also tested the model with constant term (and $k-1$ and $n-1$ dummies), which has a very low r squared (around 10%) but no collinearity problems: I would use this model if only I could separate the effect of the two reference dummies on the constant term (and I don't know how to do it).


title: Is there any optimization algorithms for these problems?

text: I have modeling the problem with the following equation: $$ \min_{X} L(X)=f(X)-\alpha g(X) + \beta k(X) $$ where $\alpha \gt 0, \beta \gt 0$, and, in my cases, I found the $f(X)-\alpha g(X) \lt 0$ and $k(X) \gt 0$, it seems this problem is not convex, how to set the parameters $\alpha$ and $\beta$? Is there any optimization algorithm for this problem? Is the Frank-Wolfe algorithm for solving these problems? Thanks!


title: How to score the options in a questionnaire and develop a scoring range in a risk profiler questionnaire?

text: I have developed a risk profiler questionnaire but don't know how to score the options and make a scoring range. Once the range is developed, I could categorize the client as conservative, moderate and aggressive. This is for a banking and broking firm . We need to analyse the risk taking capacity of the client . The questionnaire has questions which would analyse the amount of risk and the attitude towards risk . Each question has five options namely A, B , C, D ,E . The options A being the more conservative investor choices and E being the more aggressive investor choices. Kindly guide me.


title: Simulating multiple linear regression

text: I would like to simulate a multiple linear regression model using R. If I have the skewness and kurtosis for the residuals, how can I do that?


title: Given a set of data, how can I calculate the equation of the mean line that passes through the data?

text: I have a set of data that I gathered doing some tests. I plotted it and created a scatter chart. I would like to find the equation of the mean line that goes through these points. How can I do this? I remember going over this in my college stats class but have forgotten. Also it would be great if you could refer me to a tool that I could just copy and paste my data out of google docs and have it do it for me.


title: What do "endogeneity" and "exogeneity" mean substantively?

text: I understand that the basic definition of endogeneity is that $$ X'\epsilon=0 $$ is not satisfied, but what does this mean in a real world sense? I read the Wikipedia article, with the supply and demand example, trying to make sense of it, but it didn't really help. I've heard the other description of endogenous and exogenous as being within the system and being outside the system and that still doesn't make sense to me.


title: Simple mean reversion measure for binary time series

text: I am trying to define a simple measure for mean reversion in a stochastic sequence of ones and zeros, which I denote by $x_t$. Yes, a unit root test on the cumulative sum could be a viable choice, but I started off by correlating the current value of the stochastic process with the deviation of the cumulative sum from its expected trend. In more detail: let $p_1$ be the probability of $x_t=1$ (and hence $P(x_t=0)=1-p_1$), $n$ the length of the vector $x$. I then define the above mentioned deviation as ` y = cumsum(x) - (1:n)*p1 ` and my measure of mean reversion as the (negative) correlation between $x_t$ and $y_{t-1}$: ` mr = -cor(x[2:n], y[1:(n-1)]) ` Or as a function:               MeanRevert = function(x,p1){       n=length(x)       #y is the detrended random walk.       y = cumsum(x) - (1:n)*p1       mr = -cor(x[2:n], y[1:(n-1)])       return(mr)     }      Now for a purely random walk, I expected an average zero value for that correlation but there appears to be a bias:               p1 = 0.12;N=100;     MRsample=vector()      for (i in 1:500){         x= sample(c(1,0), N, replace=TRUE, p = c(p1, 1-p1))         MRsample = c(MRsample, MeanRevert(x, p1=p1))       }     hist(MRsample, xlab ="mean reversion", main =paste("p=",p1, "N=",N))      ![enter image description here](http://i.stack.imgur.com/dhsu0.png) What am I missing, what causes this bias ? Thanks! ML


title: How to use multiple datasets in order to measure the performance of a learning system?

text: I’m working on a project where I need to test a machine learning system which has a lot of hyper-parameters. Further, in order to gauge the performance of system, I’m planning to use several datasets. Firstly, I divided first dataset into three sub-sets (i.e. training, cross validation and testing). Cross validating dataset was used to find optimum values of hyper-parameters of the system. Also, I increased the training dataset size in each iteration and performance was measured using the cross validation dataset. The reason for doing this was to identify the optimum size of the training dataset. Now my question is: Do I need to do all the steps mentioned above for other dataset, or can I use the parameters (such as hyper-parameters of the learning system and training dataset size etc.) obtained above?


title: Is splitting one hurdle model in two GLM/GAM models a valid approach?

text: I came across several publications dealing with overdispersed zero-inflated count data that "simply" modelled presence absence in one model and then postive counts in a second model. This led to two models with two different outcomes. The authors stated they were using "hurdle models". In my opinion hurdle models do not work like that but "integrate" the results of both models in one step in one final output, is this right? Is the other approach now "wrong" ??


title: Does there exist some sample size threshold for choosing between the chi-squared test and the exact binomial test?

text: I'd like to use the exact binomial test to analyze some data, but the exact binomial test takes too much computation and is unnecessary when the sample size is larger because it can be approximated by the chi-squared test when the sample size is large. So what's a practical threshold for choosing between the chi-squared test and the exact binomial test? Could somebody point to some references?


title: Graphs type for likert scaled questionnaire?

text: My questionnaire is on **Likert Scale** and is divided into 7 sections say A, B, C, D, E, F, and G. In each section there are a number of questions which the respondent has to fill in accordance with the likert scale (0=No to 5=very high). To make it simpler my questionnaire looks like: > ##  **Sec: A** > >   * Q1 >   * Q2 > > > ## Sec: B > >   * Q1 >   * Q2 >   * Q3 >   * Q4 > > > And so on... ## Queries   1. I want to make **graphs** of each section (i.e. A, B, C, D, E, F, and G) in such a manner that for section "A" a combined graph is formed in which all the questions i.e. Q1, Q2, and Q3 are represented.   2. I am a civil engineer by profession therefore any guide regarding the statistics required for a research based project, will be highly appreciated.   3. Is it possible that I could upload my questionnaire here and, after seeing it any one might help me in identifying the type of statistics to be used on it?


title: How to cluster multidimensional (parametric) distributions?

text: It seems to me that it is possible to cluster unidimensional distribution in the space of their parameters, for example, we can try to find similar normal distributions in the space of $(\mu,\sigma)$ pairs. However, is there a way to cluster multinormal distributions given their means and covariance matrices?


title: Parametric spline regression after Gaussian process regression?

text: I'm trying to develop a (functional) model to predict the output of a computer simulation. I've run a small set of experiments, varying the 4 predictors $a$, $b$, $c$, and $d$ via LHS, and calculated the output, $y$. The response is expected to be nonlinear, so traditional regression seems insufficient. As a starting point, I've fit a Gaussian Process to determine as obvious relationships, see below: ![enter image description here](http://i.stack.imgur.com/kui05.png) ![enter image description here](http://i.stack.imgur.com/2XTyn.png) ![enter image description here](http://i.stack.imgur.com/pNZ0E.png) The next step is to develop the actual functional model $y \approx f(a,b,c,d)$. I've been doing a lot of reading, and it seems spline regression might be the best approach, is there a better option? If spline regression is the best approach, can anyone recommend an R package (or any other language) for this? I've tried a few, but CRAN is loaded with many "spline regression" packages, and it's difficult to identify which packages would work. Keeping in mind, I need to output the coefficients $\beta_i$ and the knots $x_i$. If my approach is completely incorrect, any help or recommendations are very welcome.


title: Relative scoring system without much effect from outliers

text: I am trying to assign scores to phone's specifications to calculate a final rating score for each phone. However, when I was scoring based on megapixels. I found that one of the phones 'A' has a 41MP camera. But all other phones have 8MP or less, which is pretty good for a phone. Right now, when I grade relatively wrt phone A (rating 5/5), all other phones get a score of 1 or 2 out of 5. However I would like it to be a little more fair, as it was just one phone with a very high value. Could you suggest me the best way to go about solving this issue where the relative grading is more fair?


title: How to justify the use of categorical variables as continuous variables in logistic regression?

text: One question again to be clarified: Can I use the variables as noted below [(3) a,b,c etc] as continuous variables in my logistic regression and if so what will be my explanation in the paper that I am writing. I have the following sets of variables:   1. A Categorical (binary) variable Ayurveda and Allopathy    2. Test variable (binary) "Spirituality is a scientific subject": Agree and Disagree    3. Then I have a number of participant perspectives/characteristics such as:      * (a) Do you believe there is life after death: 1) yes, 2) no, 3) not sure      * (b) To what extent do you consider spiritual 1) Very 2) moderate 3) slightly 4) not at all      * (c) How often would you say the experience of illness increase patients’ awareness of and focus on R/S: 1) Rarely 2) Never 3) sometimes 4) Often 5) Always 6) Not apply      * (d) etc= several more such variables with multiple choices as above Please advise.


title: How to test the equivalence of clinical parameters between two contingency tables (confusion matrix)?

text: What will be the easiest way to compare the derivations (clinical parameters) from two contingency tables? I want to test two 2x2 tables with their respective True Positive, True Negative, False Positive and False Negative entries ($H_0$: Table 1 = Table 2 _vs._ $H_1$: Table 1 ≠ Table 2). I want to compare statistically the sensitivities, specificities, negative and positive predictive values of the two tables. Sample size for the two tables is different.


title: How to compare two samples of frequencies with categorical x values where one is subset of the other

text: I am trying to compare two arrays of frequencies. The second is a subset of the first. I want to know if the second is a representative of the first. The arrays are:               x    All  Subset     a    136   38     5    127   27     b    103   23     1    102   17     6     71   11     2     27    4      I want a test that, where B behaves like A (i.e. can be viewed as A scaled), the test returns a p-value close to 1 (Null hypothesis is that for each point on x axis, percentage frequencies of B are the same of percentage frequencies of A). I tried to use Chi-squared test but, since I could have missing values, I don't know if the test validity can be compromised.. Data is missing not for technical failures but because in one subset is likely to have only few features (x values) with non-zero frequencies. Moreover since the size of the second column (the subset) is not fixed, I don't know how to scale the first column in order to obtain a valid p-value (now magnitude order is 10e-80). Thanks


title: Why use a beta distribution on the Bernoulli parameter for hierarchical logistic regression?

text: I'm currently reading Kruschke's excellent "Doing Bayesian Data Analysis" book. However, the chapter on hierarchical logistic regression (Chapter 20) is somewhat confusing. Figure 20.2 describes a hierarchical logistic regression where the Bernoulli parameter is defined as a the linear function on the coefficients transformed through a sigmoid function. This seems to be the way hierarchical logistic regression is posed in most of the examples I've seen in other sources online as well. For instance - http://polisci2.ucsd.edu/cfariss/code/SIMlogit02.bug However, when the predictors are nominal, he adds a layer in the hierarchy - the Bernoulli parameter is now drawn from a beta distribution (Figure 20.5) with parameters determined by mu and kappa, where mu is the sigmoid transformation of the linear function of coefficients, and kappa uses a gamma prior. This seems reasonable and analogous to the coin-flipping example from chapter 9, but I don't see what having nominal predictors has to do with adding a beta distribution. Why wouldn't one do this in the case of metric predictors and why was the beta distribution added for the nominal predictors? **EDIT:** Clarification on the models I'm referring to. First, a logistic regression model with metric predictors (no beta prior). This is similar to other examples of hierarchical logistic regression, such as the bugs example above: $$ y_i \sim \operatorname{Bernoulli}(\mu_i) \\\ \mu_i = \operatorname{sig}(\beta_0 + \sum_j \beta_j x_{ji} ) \\\ \beta_0 \sim N(M_0, T_0) \\\ \beta_j \sim N(M_\beta, T_\beta) \\\ $$ Then the example with nominal predictors. Here's where I don't quite understand the role of the "lower" level of the hierarchy (incorporating the logistic outcome into a beta prior for a binomial) and why it should be different than the metric example. $$ z_i \sim \operatorname{Bin}(\theta_i, N) \\\ \theta_i \sim \operatorname{Beta}(a_j, b_j) \\\ a_j = \mu_j \kappa \\\ b_j = (1- \mu_j) \kappa \\\ \kappa \sim \Gamma(S_\kappa, R_\kappa) \\\ \mu_j = \operatorname{sig}(\beta_0 + \sum_j \beta_j x_{ji} ) \\\ \beta_0 \sim N(M_0, T_0) \\\ \beta_j \sim N(0, \tau_\beta) \\\ \tau_\beta = 1/\sigma_{\beta}^2 \\\ \sigma_{\beta}^2 \sim \operatorname{folded t} (T_t, DF) $$


title: Best way to compare two alternative models

text: Let's suppose we are attempting to predict traffic accidents on some highway. Traffic Accidents would be the response variable in our data set. We have TWO competing models:   1. Two predictor variables: "Time of Day" and "Total Number of Cars Passing Particular Point Per Minute"   2. Two predictor variables: "Precipitation Level (Weather)" and "Number of Cars in FAST LANE Passing Particular Point Per Minute" The goal is to select which of these two models appears to do a better job of predicting traffic accidents. Additionally, we want to ascertain based on these two models, under what circumstances Traffic Accidents are highest and lowest. For instance, in model 1), what time of day, and what total car count per minute. Question: Which statistical software package would you use? What specific procedure would you use? Why? (Note: SAS and many other packages are available to use in our offices. SAS is preferred.)


title: how to calculate a summary value and statistical error in time series

text: I have a set of data that comes for empirical measurements over a number of days. From the beginning of the experiment to the end of it, every five minutes temperature was measured inside (Ti) and outside (To) a box, and heat flux (Q) was measured between one of the walls of the box and the outside. For every time point I use these three measures to calculate the thermal resistance, which is R = (Ti - To)/Q Obviously each R has a measurement uncertainty due to Ti, To and Q, Obviously there is correlation between measurement, because each measurement is not independent of what has been measured before. Hence my questions are: 1) what is the most appropriate way to summaries how the R value changes over time in one single value? 2) given the kind of data I describe, what is the most appropriate way of finding an error value for the value in question (1)? The time series looks like this (it has one anomalous negative peak because of some other stuff that happened to the box at the time of the peak): ![R value over time \(in seconds\)](http://i.stack.imgur.com/yLyGS.jpg)


title: Statistical Significance with large data sets

text: When I was a Ph.D. student I was trained in no uncertain terms that When we had large numbers of data results, that the number of significant results HAD to be, in and of themselves SIGNIFICANT! I am referring to the use and abuse of statistical significance in a large data set using different methodologies of predictive analytics. I would appreciate if anyone could kindly offer constructive comments with reference to this issue. Thank you


title: Compare the marginal effect for a dummy variable with a continuous variable from a probit model

text: I use probit model for one regression which consists one independent dummy variable that captures whether the firm is unionized or not. I also estimate the same probit model instead of consisting one independent dummy variable that captures whether the firm is unionized or not. This model instead have one continuous independent variable that indicates the percentage of workers who are unionized while keeping other independent variables as the previous model. I want to compare the marginal effect of a dummy variable with the marginal effect of a continuous variable, where the two variables are alternative measures of the same factor. What value of the continuous variable might allow me to compare the effect of this continuous variable with the effect of the typical union firm?


title: Sign of coefficients in linear regression vs. the sign of correlation

text: > The sign of a regression coefficient will always be the same as the sign of > the correlation coefficient between the corresponding predictor variable and > the dependent variable.   * Is the above statement always correct?   * If no, what can be the reason if the signs differ?


title: Improving a linear regression: Add predictors or change model?

text: I am trying to model a time series variable $Y_{t}$ with $4$ physical predictor variables. I used the following linear regression: $Y=\beta_{0}+\beta_{1}f_{1}(X_{1})+\beta_{2}f_{2}(X_{2})+\beta_{3}f_{3}(X_{3})+\beta_{4}f_{4}(X_{4})+\epsilon$ where $f_{i} \in \\{\log, \tan, \sin, \cos,1/.,Id\\}$. My best model gives me an adjusted determination coefficient $R^2_\rm{adjusted}$ of $0.87$. But I know that this indicator is not completely reliable. **So I would like to know if I can improve my model.** **If yes, how can I know if I have to:**   1. Change the model of regression: add some polynomials or other depedencies?   2. Look for more predictors $X_{i}$ that I may have not discovered yet? Any ideas would be nice! _**Edit: My true purpose_** $Y$ is a value which is measured each day, not at the same hour. What I really want to do is studying the variation of $Y$ with time and "cleaned" from the other external variables $X_{i}$ which are "blurring" my data. So thanks to my regression I want to find the part of $Y$ which depends on this external varying factors, and then remove this part. At the end, if my estimation $\hat{Y}$ is correct, $Y-\hat{Y}$ will correspond to my data "unblurred", and I will be free to sudy a "normalized" or "cleaned" version of $Y$ with time. I tested my model for another set of data of $Y$ measured a year ago, and it also fits the data with a $R^2_{adjusted}$ around $0.8$. And I can see graphically that the model is not completely insane, but not perfect either.


title: Are pairwise Wilcoxon tests a valid non-parametric alternative to Tukey's HSD test?

text: I am trying to identify significant differences between groups from normal and non-normal distributions. When the distribution is normal the workflow is pretty straightforward: ANOVA and after that a Tukey's HSD test in order to identify the groups which are responsible of the significant variation. When the distribution is non-normal I apply a Kruskal-Wallis test and after that a pairwise Wilcoxon test (R function `pairwise.wilcox.test()`) to identify the groups. I know that Wilcoxon test is the non-parametric version of the _t_ test, and is not corresponding to the Tukey's HSD test, but I want to know if what I am doing is wrong and which is the main difference between the outcomes of a _t_ test and a Tukey's test since they both give differences between groups.


title: Sample size for a non-normal distribution

text: I'm quite new in this field. I hope my question makes sense. I have a database that stores information for around 10.000.000 projects. Each project has several features (let's call them X) like number of lines of code, number of people, etc. In addition, there are some extra features (let's call them Y) that I cannot calculate directly on the 10.000.000 projects due to time and storage limitations. Therefore, I would like to take a representative sample of these projects and perform a qualitative analysis to calculate the features Y. The frequency distribution for each of the X features seems to follow an exponential distribution (I run the one-sample Kolmogorov-Smirnov test and I got a p-value < 0.001). In addition, I generated some random samples (with sizes of 100, 500, 3000 projects) to calculate the features Y, and even those features seem to follow a negative exponential distribution. I would like to know if there is a formula to get the minimum sample size to say that my sample is representative enough. Can I use the concept of confidence level (95%) and interval of confidence (1.96) to calculate the minimum sample size even if my distribution is not normal? thanks in advance


title: Representative point of a cluster with L1 distance

text: The representative point of a cluster or cluster center for the k-means algorithm is the component-wise mean of the points in its cluster. The mean is chosen because it helps to minimize the **within cluster variances** (which is to say that it is minimizing within cluster squared Euclidean distance, since its the same). Now, when I want to use a different distance metric, like L1 (Manhattan distance), I would be minimizing the **within cluster absolute deviations** (if I choose to). I was looking at what the centers calculation would be like, because the mean won't reduce the new objective function, and I found this on the Matlab documentation of their `kmeans` function for the distance parameter > `cityblock` \- Sum of absolute differences, i.e., the L1 distance. Each > centroid is the component-wise median of the points in that cluster. Now, how did they come up with the component-wise median? How does a median minimize absolute deviations?


title: Pooled vs Unpooled variance for AB Testing

text: I'm trying to A/B test two proportions and calculate the Z-score. My hypothesis is that the two groups are the same. Most of the formulae show the use of unpooled variance. Should pooled variance be used instead? Thanks https://developer.amazon.com/sdk/ab-testing/reference/ab-math.html http://visualwebsiteoptimizer.com/split-testing-blog/ab-testing-significance- calculator-spreadsheet-in-excel/


title: Reducing variables and correlated variables in a dataset

text: I am attempting to produce models on a dataset with approximately 60 variables and 800 subjects. Other datasets are also being considered for use. The datasets were not originally constructed for this type of analysis. The dataset also comprises a mixture of variables. Some are continuous and some discrete. Also some of these are correlated with one another as can be seen on a correlation matrix and also just by eye/common sense. Currently clustering variables has been used with the aim to remove any correlated components being used to generate models to avoid over-fitting and reduce the number of variables. By choosing the most representative variable from each cluster however I am unclear if this is being done correctly and not confident this is correct and will result in incorrect models. Previously what has been done is to just load all the variables into the platform and produce a cluster summary. The most representative member variable is picked out of each cluster and models made. In many clusters there is a variable present that shows a low R squared with the cluster it is in, yet it has been included. The variable also does not show correlation to the other group members on reference to the correlation matrix. It's almost as if it has been forced in where it will fit best rather than left on its own. (Sometimes a cluster does include 1 variable.) Also, sometimes two separate clusters appear to contain variables which are correlated. Thus I feel at times correlated variables are being used for models and at the same time uncorrelated variables are being missed. I am very new to this level of stats and I am struggling with wrapping my head around eigenvalues etc. It was implied that this is an aspect of PCA but I have found that this is not really the case. Information on the web tends to get into some quite complex math on the subject. Here is my understanding, it may well be off. I understand that PCA is a transformation of the dataset so that the first component explains the most variation in the variables, the second PC is orthogonal to this (uncorrelated) and that there is one PC for each variable. Each PCA is given a total eigenvalue that = number of variables (due to the transformation of the original axis where the variance of the data is made to =1 and mean =0). The PC that explains most of the variation in the data takes the greatest share of this followed by the second etc. This is due to variable weighting coefficients where each variable is given a weighting in respect to its influence on that PC. Eigenvectors are where I fall apart!! Clustering involves looking at all the variables and performing iterative splits. A split involves a PCA analysis on the cluster members? then.. from a guide.. 'the cluster with the largest second eigenvalue is chosen to be split into two new clusters' Again I fall apart due to not really understanding eigenvalues. I think the tools I have been using have been suggested to perform a dump and click solution which I can't really believe..but I cannot argue the case due to limited understanding. Long post, but I'm getting a bit lost!


title: How to deal with nominal variable with too many levels?

text: currently I'm trying to model a response variable y, and I have zip code as my independent variable, my model is logistic regression. When it comes to nominal variable, the text book method is to create k-1 dummy variable (assuming the nominal variable have k different levels), but zip code's k is too big, I can't create that amount of dummy variables, is there any other ways to deal with this? Or more generally, how to deal with nominal variables with too many levels (k>=100)?


title: Should I include multiple mediating variables in the same analysis?

text: If the answer is "it depends", what does it depend on? I have three mediating variables, one predictor variable and one outcome variable. I have observed that the results I get when I run three separate analyses are different from those I get when I put them all in the same analysis. Why is that? I'm using Andrew Hayes's PROCESS syntax for SPSS.


title: Understanding standard errors in logistic regression

text: I am having trouble understanding the meaning of the standard errors in my thesis analysis and whether they indicate that my data (and the estimates) are not good enough. I am performing an analysis with Stata, on immigrant-native gap in school performance (dependent variable = good / bad results) controlling for a variety of regressors. I used both logit and OLS and I adjusted for cluster at the school level. The regressors which are giving me trouble are some interaction terms between a dummy for country of origin and a dummy for having foreign friends (I included both base-variables in the model as well). In the logit estimation, more than one of the country*friend variables have a SE greater than 1 (up to 1.80 or so), and some of them are significant as well. This does not happen with the OLS. I am really confused on how to interpret this. I have always understood that high standard errors are not really a good sign, because it means that your data are too spread out. But still (some of) the coefficients are significant, which works perfect for me because it is the result I was looking for. Can I just ignore the SE? Or does it raise a red flag regarding my results? I usually just ignore the SE in regressions (I know, it is not really what one should do) but I can't recall any other example with such huge SE values.


title: Evaluating productivity metrics in an A/B test system

text: I'm beginning work on an A/B testing system for UX changes to a tool my team builds. Our users are doing simple, repeatable tasks assigned by another team's system. My team handles the UI/UX part. Users typically have a steady stream of work assigned to them, and they can see their current task, and complete it. Tasks typically take under a minute. Each task has an expected time to complete it and the users are rated based on their overall "actual / expected" for a day's work. My question is this: Given that a random subset of users are consistently seeing a slightly different treatment (the deterministic 'random' selection is a solve problem), and I can get a steady stream of [User, expected time, actual time, treatment] for each complete task, how can I evaluate the overall productivity changes for users in a statistically valid way? To give some insight into why this is tricky, here are some things we've been worried about:   1. Let's say users under Treatment A are 10% more efficient. As a result, they will produce _more data points_. With the wrong analysis, more data points might skew the results, and we want exact figures, with high confidence. We want to know how much the change helps the average user.   2. The data is very noisy, but seems to average out over an hour/day/shift for each user. One task may take a fraction of the expected time, or a magnitude more than expected. This is in part because our ability to estimate completion time is in need of improvement (but on average, it's pretty accurate). We worry that our high variability in estimated completion time will increase our uncertainty in results analysis.  Thanks very much.


title: help with AUC for PR curve when data has tied and skewed values

text: I am wondering if there are methods available to calculate AUC for Precision Recall curves when the predicted scores/probs/beliefs(whatever you want to call it) has tied values and could be skewed too. I think data being skewed is a subcase of having tied values. I have seen the following paper which talks about not doing linear interpolation(section 4) but doesn't go into the possible issues with tied scores. The Relationship Between Precision-Recall and ROC Curves For Example               ID      BELIEF     P1      0.1     P2      0.1     P3      0.1     P4      0.1     P5      0.5     P6      0.5     P7      0.5     P8      0.8     P9      0.8     P10     0.8      -Abhi


title: Feature importance and interpretation of alternating decision trees

text: Is there a way of calculating feature importance in alternating decision trees? What if I've already trained an alternating decision tree and want to calculate feature importance in terms of classifying a specific instance. In my understanding, the overall prediction is calculated by summing the contributions of all the prediction nodes through which the instance passes. Would a reasonable measure of feature importance for that particular instance just be the sum of nodes involving that feature?


title: Problems with SEM: Non-positive definite matrix

text: Link to picture: https://www.dropbox.com/s/bxzxy0y37npqc4o/SEM%20problems.PNG Need a bit of help determining what I should do to fix this problem: The SEM wont run because it believes it's non-positive definite. Meaning that there may be singularity. However these items are very different as i have items measuring autonomy, job feedback, coworker support, positive mood and employee engagement. Can anyone suggest a solution to this problem?


title: Finding appropriate distribution that fit to for a frequency distribution of a variable

text: I am trying to fit distribution to my dataset. For example I have data whose frequency distribution look like the following. ![enter image description here](http://i.stack.imgur.com/rvF4E.jpg) The following is frequency table, after divided the data in 10 bins.                 xvar lvar yvar     1     1   y1   25     2     2   y1   22     3     3   y1   16     4     4   y1   15     5     5   y1   10     6     6   y1    2     7     7   y1   NA     8     8   y1   NA     9     9   y1   NA     10   10   y1   NA     11    1   y2    2     12    2   y2    5     13    3   y2    8     14    4   y2   10     15    5   y2   15     16    6   y2    9     17    7   y2    6     18    8   y2    3     19    9   y2    2     20   10   y2    1     21    1   y3   NA     22    2   y3   NA     23    3   y3   NA     24    4   y3    6     25    5   y3   12     26    6   y3   14     27    7   y3   18     28    8   y3    6     29    9   y3    2     30   10   y3    1      `yvar` is frequency in each bin. I would like to find a corresponding distribution that best fits this trend (may be beta distribution or normal distribution?). How can I achieve this? **Edits** : (1) I do not need perfect fit but close fit and simple: A commonly used distribution is good enough (2) The `xvar` is a continuous variable (like blood pressure) or discrete scale (for example 1:10 rating scale) (3) I believe that most of it can be described by normal distribution but I have some of them skewed like `y1`. I wonder such values can be described by beta distribution, but do not know how to find appropriate value (beside just guess)? The solution does not necessarily need to be in `R`, but following is `R` code for the above data set and plot.               y1 = c(25, 22,16, 15, 10,   2, NA, NA, NA, NA)     y2 = c(2, 5,8, 10, 15, 9, 6, 3, 2, 1)      y3 = c(NA, NA,NA, 6, 12, 14, 18, 6, 2, 1)            myd <- data.frame (xvar = 1:10, lvar = rep(c("y1", "y2", "y3"), each = 10),  yvar = c(y1, y2, y3))      library(ggplot2)     ggplot(myd, aes(x = xvar, y = yvar, fill = lvar, colour = lvar)) +  geom_point () + geom_line ()  + theme_bw()


title: Longitudinal modelling of ordinal censored data

text: I have n subjects observed over some time $t{\scriptscriptstyle n}<T$ . At several time points during $T$, I measure in each subject the value of an ordinal variable $V$ (which can take 5 values). Not every subject has measurements in every time point and there is censoring of the subjects during T. I would like to model the probability of a random subject being in each of the categories of V over time, allowing for a random subject effect. What would be the most appropriate methods to perform this? Is it possible in R and could someone pinpoint me to any information sources about it?


title: Series dimensionality reduction for classification Input

text: I am looking to construct a predictive model where the outcome variable is binary and the input is time series. To make it more concrete, the model will predict if a customer churns (left the company; coded as 1 or 0) based on the amount they spent with the company in the prior 60 days. So, the data is one customer per row and the columns are an outcome factor (1 or 0) and 60 additional columns for the amount spent the in time t-1, t-2....t-60. Here is some example data:               #create the data a series of length 60 and a class ID     sc <- read.table("http://kdd.ics.uci.edu/databases/synthetic_control/synthetic_control.data", header=F, sep="")          #binary class lable     classId <- as.factor(c(rep(0,300), rep(1,300)))     newSc <- data.frame(cbind(classId, sc))     newSc$ID<-seq(1,600,1)      The actual model may have many of these series for each customer, so I need to reduce the dimensionality of the data for the series, e.g. instead of using 60 values, I need to reduce this down to a handful. Of course, I can use the mean, min, max etc of the series but I have been reading about using Discrete Fourier Transform. **Questions:**   1. Is the DFFT in R a proper method to use for my purpose? Any information on how it works would be appreciated.   2. Assuming this R function is correct, how do you extract just the most meaningful coefficients to achieve dimensionality reduction? ADD: There seems to be a consensus that using DFFT for dimension reduction is not a wise choice, but it seems that in data mining, this function, DWT and SVD are all commonly used: Time Series Mining starting on page 20.


title: To exclude or include the intercept in GLM model

text: When is is appropropriate to include or exclude the intercept from a regression model? SPSS provides this option in the GLM menu. I am assessing group difference (gender) on tasks performed (X) and income earned (Y). In this case, if no task is performed, then no income is earned. The SPSS guidelines states here: "Include intercept in model. The intercept is usually included in the model. If you can assume that the data pass through the origin, you can exclude the intercept." My stats knowledge is a bit rusty now but I remember my teacher saying that you should never exclude the intercept!


title: Effect size as the hypothesis for significance testing

text: Today, at the Cross Validated Journal Club (why weren't you there?), @mbq asked: > Do you think we (modern data scientists) know what significance means? And > how it relates to our confidence in our results? @Michelle replied as some (including me) usually do: > I'm finding the concept of significance (based on p-values) less and less > helpful as I continue in my career. For example, I can be using extremely > large datasets so everything is statistically significant (p<.01) This is probably a stupid question, but isn't the problem the hypothesis being tested? If you test the null hypothesis "A is equal to B" then you know the answer is "No". Bigger data sets will only bring you closer to this inevitably true conclusion. I believe it was Deming who once gave an example with the hypothesis "the number of hairs on the right side of a lamb is equal to the number of hairs on its left side." Well, of course it isn't. A better hypothesis would be "A does not differ from B by more than so much." Or, in the lamb example, "the number of hairs on the sides of a lamb does not differ by more than X%". Does this make sense?


title: Multiple testing correction

text: I use `glm.nb` from `MASS` package over many observations. I then extract the p-values for the interaction term. Due to multiple testing, I use BH method for correction for multiple testing, with FDR of 0.05. Sometimes, when I apply the negative binomial model to fit to the data, the model doesn't converge and ends up with a warning "Iteration limit reached" and has a huge theta value normally (> 10000). However, most of these times, the p-values obtained are quite small (< 10e-5). Of course these small values are wrong, as there was no convergence achieved. Now, I can find out these observations using the $th.warn attribute which stores the warning message. My question is, before performing multiple testing, should I remove these non- fit observations and therefore wrong p-values (as they might affect the BH process of correcting for FDR) or because I have already performed the test, removing these entries after finding the warning message would bias the FDR correction of the p-values and therefore I should just continue with FDR correction without removing these obviously wrong but small p-values??? Thank you!


title: Parameters to change for different kernels for SVM

text: I am carrying out SVM and was interested in knowing the parameters that could be varied for each kernel. I am using 3 kernels: RBF, linear and polynominal. These are the parameters that i think can be changed:               For Linear : C     For RBF: C and gamma     For Polynominal : C and gamma      Would like to hear some feedback whether i am wrong or correct. If there are any other parameters, do tell me... Thanks.


title: Breakpoint for bivariate data

text: The breakpoint(s) estimation approach implemented in the strucchange package (Zeilei & al) seems to work very well (based on my little experience with this package on real case studies). Is there an existing similar approach for multivariate regression models ? Otherwise, would it be easy to generalize the strucchange method ? (Unfortunately I do not have time to study this question myself). ## Update Here is an animated example (thanks to Yihui Xie's animation package). I have a bivariate dataset, each pair $(x,y)$ has been recorded at some date:                        x        y       date     1 3.690131 6.797999 25/01/2012     2 3.552278 7.055363 07/02/2012     3 3.623821 6.745984 10/02/2012     4 3.600735 6.847450 10/02/2012     5 3.726609 6.894321 14/02/2012     6 3.578204 6.823344 17/02/2012     ...      Sarting $n_1=40$, I estimate the shape of the isodensity ellipses of the first $n_1$ observations as well as the shape of the isodensity ellipses of the remaining $n_2=n-n_1$ observations. These are simulated data. For some $n_0>40$ I have simulated Gaussian bivariate i.i.d pairs $(x_i,y_i)$, $i=1, \ldots, n_0$, and then I have simulated Gaussian bivariate i.i.d pairs $(x_i,y_i)$, $i=n_0+1, \ldots, n$ with different parameters. The question is: find $n_0$, which is here a breakpoint for the model `lm(cbind(x,y)~1)` in R syntax. ![enter image description here](http://i.stack.imgur.com/4HSSp.gif)


title: How to compare two regression slopes for one predictor on two different outcomes?

text: I need to compare two regression slopes where:               y ~ a + b1 x     y ~ a + b2 x      How can I compare b1 and b2? Or in the language of my specific example in rodents, I want to compare               antero-posterior diameter ~  a + b1 * humeral length        de naso-occipital length  ~  a + b2 * humeral length


title: sequence of regularization values in cross-validation

text: Does anybody know how cv.glmnet (R's package glmnet) or LassoCV (scikit-learn) chooses a sequence of regularization constants, which it uses in cross- validation? Are there any recommendations on how to pick a range of regularization values in lasso, ridge regression and elastic net? Thank you very much!


title: Quick question regarding simple algebra in Bayesian statistics

text: In this paper you can see equation (11), the _probability density_ , written as: $$S_n(\xi) = \frac{1}{N}\frac{1}{2\pi\sigma_r(n)\sigma_g(n)} \sum\limits_{i=1}^{N}\, \exp\left[-\frac{1}{2}\left(\frac{x_i^2}{\sigma_r(n)}+\frac{y_i^2}{\sigma_g(n)} \right) \right]$$ never mind what each parameter represents, my question regards the next two steps, where it first defines the likelihood (eq. 12) as: $L=\prod\limits_{n=1}^{M}\, S_n(\xi)$ and finally it applies a logarithm to that last equation to obtain (eq 13): $$\ln \,L = -\sum\limits_{n=1}^{M} \ln\left\\{N2\pi\sigma_r(n)\sigma_g(n)\,-\,\sum\limits_{i}\exp\left[ -\frac{1}{2}\left(\frac{x_i^2}{\sigma_r(n)}+\frac{y_i^2}{\sigma_g(n)}\right) \right] \right\\}$$ Now, I'm either missing something **really** obvious here or that last equation is **wrong** and should be: $$\ln L = -\sum\limits_{n=1}^{M} \ln\left\\{\frac{N2\pi\sigma_r(n)\sigma_g(n)}{\sum\limits_{i}\exp\left[ -\frac{1}{2}\left(\frac{x_i^2}{\sigma_r(n)}+\frac{y_i^2}{\sigma_g(n)}\right) \right]} \right\\}$$ or equivalently: $$\ln L = -\sum\limits_{n=1}^{M} \left\\{\ln(N2\pi\sigma_r(n)\sigma_g(n))\,-\,\ln\left(\sum\limits_{i}\exp\left[ -\frac{1}{2}\left(\frac{x_i^2}{\sigma_r(n)}+\frac{y_i^2}{\sigma_g(n)}\right) \right] \right) \right\\}$$ Is equation (13) in that paper wrong or is my algebra terrible?


title: Determining product sales from incomplete sales information

text: Apologies if: a) This is a bit convoluted; or, b) this would be more appropriate on a different Stack site (though I don't think the latter is the case). My company has two categories of products we sell, call them A and B. Within each category, there are 9 products and 7 products, respectively. So, A1:A9 and B1:B7. These products are sold through two types of terminals: a cash-register style terminal at retailers (T), and in vending machine-style machines (V). The V terminal records all transactional details- I know the amount of money put in, the exact mix of product sold from both categories, cashes (or returns) of all products, and any vouchers printed out for unspent money. From the moment money is entered into the machine, I know every selection made by the customer until all the money for that session is used up. The T terminal is a different story. It tracks none of the money put in, and doesn't track sales of Category B products. It only tracks returns of all products and sales of products A1:A9. The only way we know sales of Category B products from Terminal T is by looking at weekly, monthly, or yearly settlements (where the retailer actually paid for the product after it was sold). To recap: I know: Terminal V: Cash spent, products sold from both categories, and products returned for credit for each customer using the terminal. Terminal T: Products sold from Category A and products returned for credit from both categories, for each customer using the terminal. Weekly/Monthly/Yearly: The amount of Products B1:B7 sold by the retailer at all available terminals. Customer Behavior at both Terminals: We know that a high percentage of people returning product spend the money from the return on new product, rather than pocketing the cash (We get this info from Terminal V and assume it also applies to Terminal T). What I'd like to model from this information: How much customers are spending per session on Category B at Terminal T. With what I know, what would be a good way to mathematically maximize a model for determining player behavior when it comes to the purchase of Category B at Terminal T? Please let me know what clarification would be helpful. I know this is kind of a mess.


title: Block bootstrap for dependent data with unequal sampling intervals?

text: I have data from a natural archive (lake sediment). For various reasons it is usually impossible to sample the archive equally in time, and we end up with a time series where essentially we have observations at ~~random~~ [*] time points. (In other words we are not sampling a regular process and missing some of the regular observations - we have a continuous process [sediment accumulates constantly over time] and we observe something from that process essentially ~~randomly~~ in time.) For normal time series data one might employ a block bootstrap when applying the bootstrap to problems of confidence interval generation etc. Does the block bootstrap depend on equal sampling intervals or would it be appropriate to apply it to the unevenly-sampled data I describe? [*] Apologies, the use of "random" here is wrong. We sample the mud evenly in terms of depth (or distance from the mud-water interface, larger distances == older samples). However, the accumulation rate of the sediment is rarely constant, and even if it was there is the issue of compaction of the sediments, hence each, say 0.5cm, slice of mud represents a slightly different amount of time. Then, in general, we might only analyse at best every other slice of mud, and usually at coarser resolutions. Hence the unequal sampling _in time_ of our resulting observations.


title: Resampling and multiple comparisons

text: I hope someone can help me with the following problem (please have in mind that I am not a statistician). I run an experiment with 4 treatments where I measured the variable “P”, which can be either zero or a positive value. From each treatment I calculate an index “E” (the average of all “P” values larger than zero multiplied by the ratio of “P” values larger than zero to the total number of replicates). Therefore, for each treatment I get only one value of “E”. To present the results from this experiment I used the following resampling procedure:   * I created 1000 different sets by means of resampling with replacement from the data for each treatment   * For each resampled set, I calculated index “E” again, hence having 1000 “E” values for each experiment   * Then I calculated the average and the 95% confidence intervals for these 1000 “E” values My question is: How can I make multiple comparisons (among treatments) using the resampled “E” values? My first idea was to test the resampled data using multiple t-tests and correct the p-value using Bonferroni correction, but I am not sure it is ok. Thanks!


title: A good alternative to data binning?

text: I read many times that data binning of continuous variables is a very bad idea. For instance, let's take something like heart rate and let's define the following 2 bins: (125 - 135), (136 - 145) Let's say that (136 - 145) corresponds to a hard effort. If your exercise session causes your heart rate to stay at 135 consistently, data binning will reveal that you spent no time exercising hard, while you were 1 beat per minute away from that bin during the whole time. Obviously, this is an exaggerated example to illustrate the point. I was thinking to weigh each second spent at a given heart rate based on its distance from the centers of the bins it falls between. For instance, in the example above, 140.5 would be 100% in the second bin, 135.5 would be 50% in the first bin and 50% in the second bin, 130 would be 100% in the first bin. Does that sound like nonsense or a reasonable solution? What would be a better way?


title: Variable selection in large datasets

text: I'm looking for an overview of some methods of variable selection. I use datasets with around 6000 variables (the level of missing values is satisfying i.e. there are no variables with 100% missing values). There is no problem of insufficient degrees of freedom in my datasets. In order to use some data mining techniques (such as regression, gradient boosting, etc) I use some methods of removing "weak" variables. I am aware of stepwise methods for regression, but I'm looking for use with broader range of methods. Does anyone know of any scientific papers or other articles about this subject?


title: Ancova without intercept in R

text: Here is an exemple of fitting an ANCOVA model without intercept in R:               > k <- 40     > x <- 1:k     > sigma <- 30     > y1 <- rnorm(k, 1*x, sigma)     > y2 <- rnorm(k, 2*x, sigma)     > y3 <- rnorm(k, 3*x, sigma)     > x <- rep(x, 3)     > y <- c(y1,y2,y3)     > group <- gl(3,k)     > summary(lm(y~0+x+group:x))          Call:     lm(formula = y ~ 0 + x + group:x)          Residuals:         Min      1Q  Median      3Q     Max      -71.110 -20.389  -2.496  15.913  76.937           Coefficients: (1 not defined because of singularities)              Estimate Std. Error t value Pr(>|t|)         x          3.1933     0.1989  16.054  < 2e-16 ***     x:group1  -2.0444     0.2813  -7.268 4.44e-11 ***     x:group2  -1.2994     0.2813  -4.619 9.97e-06 ***     x:group3       NA         NA      NA       NA         ---     Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1           Residual standard error: 29.6 on 117 degrees of freedom     Multiple R-squared: 0.7654,     Adjusted R-squared: 0.7594      F-statistic: 127.3 on 3 and 117 DF,  p-value: < 2.2e-16       The parameter "x:group3" is not estimated. Of course this parameter actually is the parameter "x". The estimates seem to be correct, but are the p-values correct ? And is there another to fit this model without obtaining this problem ?


title: 2 correlated variables , 2 conditions and 3 groups: what analysis?

text: In a design with two correlated dependent variable measures with similar range and variance, two within subject conditions and three diferent participant groups each tested at a diferent time point. Which would be the most apropriate statistical analysis? I thought about a repeated measures anova, but I am not sure whether to run separate analysises or how to analyse the two correlated variables? I am interested in possible interactions between factors and also whether the factors predominantly effect just one or both of the variables. I am using SPSS 18, for the analysis.


title: Loss functions, decision theory for hyperparameters, or estimating the variance of an unknown prior

text: **How I got here** I was originally interested in **point estimation** and the Bayes risk of some distribution $\pi$ $$ r(\pi) = \mathbb{E}_{\pi(x)}[ \mathbb{E}_{\Pr(y|e,x)}[L(x,\hat x(y|e))]], $$ where $L$ is some loss function and $\hat x$ is the posterior mean (here: Bayes estimator). Here $e$ is some experimental design parameter. If $L(x,\hat x) = |x-\hat x|^2$ (squared error), then the Bayes risk is also identical to the expected posterior variance. So to optimally design experiments, I can use **negative expected variance** as a utility function to maximize. Now suppose I iterate this procedure (sequential experimental design). Although it is not necessarily the globally optimal solution, minimizing expected loss at each step is locally optimal (the greedy solution) and works quite well when compared to any non-optimized experiment designs. **The problem** In practice I need to apply these sequentially designed experiments to members of an **ensemble each with a different set of parameters $x$** whereas originally I assumed the parameters $x$ are chosen from the prior $\pi$ and fixed throughout. Now, clearly point estimation is not the appropriate interpretation of this problem. However, assuming the support of true distribution is contained in that of the prior I chose, my posterior should converge. But it doesn't! And the reason is because I am adaptively minimizing the expected variance of the posterior assuming the parameters are going to be the same at the experiment. **The solution** There are two types of solutions I can foresee being given (but I of course welcome a different type):   1. You are using the wrong loss function. You need to introduce a hyperprior on some hyperparameters, define a loss on them and try to design experiments to minimize the Bayes risk. And to this I ask: are there any treatments in the literature dealing with this decision theoretic framework for hyperparameters, or more specifically is there any natural loss functions on hyperparameters which reflect the error in estimating the variance of an unknown hyperprior?   2. Your decision theoretic framework is not really suited to your problem. You should look at _these_ papers or search Google for _this_ jargon.


title: Support Vector Machines (SVM) maximum margin hyperplane use

text: I have been approaching to Support Vector Machines in this period, and I need a clarification about the use of the maximum margin hyperplane. From a labeled point training set, we train a model that will correctly classify the future test set points. To do this, we search for a best hyperplane that divide the training set points into two classes (True, False). So we use kernel and etc to get this goal. But, when we've gotten our best hyperplane, what do we do with it? Do we just use it to classify new data points, by verifying which ones are "on this side" and which ones "on the other side"? I've understood how to get the famous best hyperplane, but I need some clarification on what to do next with it... Thanks!


title: Method to output statistical score when validating using bootstrap method on ordinal logistic regression model

text: It is about sas programming. Thank you in advance. I build an ordinal logistic regression model on my data.               Level  Meaning       3     Alert       2    Unusual Event       1    Non Emergency      As you can see, it is my ordinal response variable. 3 means large loss and Alert situation, etc. I build a very simple model, the graph below shows the result:                proc logistic data = sasuser.Pacific_West_ outmodel=sasuser.Pacifici_West_model;;          model Numerical_Emergency(event = '1') = AGE;      run;      ![enter image description here](http://i.stack.imgur.com/xBihb.jpg) Everything looks fine. Now I would like to apply the BOOTSTRAP method to validate this model. What thing I would like to see is the test c value as ROC value. As you may see in the graph above, the c (ROC) value is 0.640. The basic idea of bootstrap is to generate new data set by random picking process from the data we used to train the model as we all know. I prefer to generate 1000 copies of data with same size of original data for building the model in graph above. THE PROBLEM: In the SAS Score statement, I just found outroc method which does not support the ordinal regression. Is there a simple procedure to help me to output the statistics about fitting of trained model on the validation dataset build by myself? Thank you so much for your time. Have a nice day! Wenhao SHE


title: How to compute expectation for poisson distributed variable in context of binomial?

text: I'm working on a problem that needs to be solved using EM algorithm. In doing that, I have to evaluate an expectation that I actually have no idea how to. Consider: $Y$ as a a fixed observed integer and consider $Z$ as a random variable which has a Poisson distribution with rate $\theta$: $$ Z \sim \operatorname{Pois}(\theta)$$ $$ Z|Y \sim \operatorname{Binom}(n = y, p)$$ How can I compute the expectation below with respect to $Z$? $$E\left(\frac{y!}{Z! (y - Z)!}\right)$$


title: Using Adaboost with SVM for classification

text: I know that Adaboost tries to generate a strong classifier using a linear combination of a set of weak classifiers. However, I've read some papers suggesting Adaboost and SVMs work in harmony (even though SVM is a strong classifier) in **certain conditions and cases**. I'm not able to grasp from an architectural and programming perspective how they work in combination. I've read many papers (maybe the wrong ones) which didn't explain clearly how they work together. Can someone throw some light on how they work in a combination for effective classification? Pointers to some papers / articles / journals would also be appreciated.


title: regarding conditional independence and its graphical representation

text: When studying covariance selection, I once read the following example. With respect to the following model: ![enter image description here](http://i.stack.imgur.com/EftcN.png) Its covariance and inverse covariance matrix are given as follows, ![enter image description here](http://i.stack.imgur.com/NJ4gb.png) I do not understand why the independence of x and Y are decided by inverse covariance here? What are the mathematical logic underlying this relationship? Also, the left graph in the following figure is claimed to capture the independence relationship between X and Y, why? ![enter image description here](http://i.stack.imgur.com/kOtw8.png)


title: error of the mean in presence of background

text: Suppose I have a normal distribution, $N[\mu,\sigma]$, and I have a sample of size $n$. It is well know that the error (std deviation) of the mean is $\sigma/\sqrt{n}$. Now suppose that my distribution is a mixture of signal (normal distributed) and a flat background $$\text{pdf}[x|\mu, \sigma, s, b] = \frac{sN[x|\mu,\sigma] + b U}{s+b} $$ suppose that the uniform distribution is much wider than the normal but in the signal window the background is important. Question: what is the error of $\mu$ now? ($s$ and $b$ are random variables, they can be assumed to be Poisson distributed)


title: Regression estimate of a non-negative variable

text: I have to estimate linear weight $\beta$ for regression $Y \sim \mathbf{X}$, where $Y$ are non-negative samples. If I perform vanilla regression (lets assume ridge regression) it will find $\beta$ such that most of estimated $\hat{Y}$ are non-negative. But some of them will be negative. But I need them to be invariably non-negative. In Bayesian framework, it implies I need a prior $\beta$ that is non-negative (ridge assume normal priors). Is there simple way to modify ridge regression to obtain this? Or putting non-negative priors and estimating via MCMC is only solution to it. For further description of problem (why I need estimated $\hat{Y}$ to be positive). My full model is a hierarchical model where there is a variable $Z$ that is a linear model of latent variable $Y$ i.e $Z \sim N(\Theta^T \mathbf{Y}, \sigma^2)$. Observable variable are $Z$ and $X$. I need estimates of $\beta$ and $\theta$. I am currently estimating $\beta$ given Y. Given $\hat{Y}$ I am estimating $\theta$. I am iterating between this 2 steps. If $\hat{Y}$ is negative the problem become unidentifiable.


title: Challenging Likelihood Ratio Test

text: Derive Likelihood Ratio Test of size $\alpha$. H$_0$: $\theta=\theta_0$ H$_1$:$\theta \neq \theta_0$ \begin{equation} {f}(x,\theta , c) = \theta(x-c)^{\theta-1} \ {c<x<c+1} \end{equation} I am having trouble with the mle for c. We have to use order stats and I am not sure which is the correct one to use


title: Display Correct Poisson Calculation

text: ![enter image description here](http://i.stack.imgur.com/Huvtw.png) I'm having difficulty calculating the desired result (i.e. `0.0723`) on my calculator. Here is what I have tried:               E^1,8 = 6,049647464          (1,8^4) / 4 = 2,6244          6,049647464 x 2,6244 = 15,8766948


title: Entropy of Inverse-Wishart distribution

text: What is entropy of Inverse-Wishart distribution? I need just a reference, but derivation (e.g. using inverse property) would be interesting too.


title: 1-d extension of 2-d variance?

text: A random variable X={X1..Xm} and Var(X) is the mXm variance-covariance matrix. Is there an accepted 1-d statistic like variance that may be extracted from the sample variance of data representing a 2-d random variable? (explanation of the 2d equivalent of a sample variance) How to find variance between multidimensional points?


title: Experimental design type?

text: A pencil depth tester is used for various pieces of paper. Two different pencils are used. The class measures how deep the pencils go into the various pieces of paper. Paper Pencil Depth               1      2        5     2      1        3     3      1        2     4      2        6     5      2        7     6      1        8     7      2        2     8      1        9      What type of experimental design is used?


title: How should I approach modeling these subjective probability estimates?

text: In my data about 1000 people have made estimates of the probability of 100 unique events. On average people forecast on about 50 events, but some forecast on all of events and some on only a few. About 500 of these people are working in teams of 15. They see each others forecasts and communicate with each other, but each individual still offers their own estimates). One question I'd like to answer is how much of the variance in these forecasts is explained by which team people were (randomly) assigned to. But I also suspect that the share of the variance explained by team will vary by event. So I'd like to examine that variation as well. I have a bit of experience with mixed models with lmer, but its not immediately obvious to me how to approach this.


title: How to fit a model to self-reported number of friend interactions over a 20 day period?

text: I am a novice in statistics so please correct me if I am doing something fundamentally wrong. After wrestling for a long time with R in trying to fit my data to a good distribution, I figured out that it fits the Cauchy distribution with the following parameters:                  location      scale         37.029894   18.678936       ( 3.405665) ( 2.779136)      The data was from a survey where 100 people were asked how many friends they talked to over a period of 20 days and I am trying to see if it fits a known distribution. I generated the QQ-plot with the reference line and it looks like the image given below. From what I have been reading on the web, if the points fall close to the reference line then it is a good evidence that the data comes from this distribution. So, is this a good evidence to say that the distribution is Cauchy or do I need to run any more tests? If so, can someone tell me the physical interpretation of this result? I mean, I read that if the data falls into a Cauchy distribution, then it will not have a mean and standard deviation but can someone help me understand this in plain English? If it does not have a mean then from what I understand, I cannot sample from this distribution. What is one supposed to infer about the population based on this result? Or should I be looking at other models? **UPDATE: What am I trying to achieve?** I am trying to evaluate how much time it takes for some arbitrary piece of information to propagate for a population of size X. As this depends on the communication patterns of people, what I was trying to do was to build a model that could use the information from the 100 people I surveyed to give me patterns for the X number where X could be 500 or 1000. QQ-Plot ![alt text](http://i.stack.imgur.com/HvmjM.png) Density Distribution of my data ![alt text](http://i.stack.imgur.com/Ob7hS.png) Cauchy Distribution ![alt text](http://i.stack.imgur.com/qGHmd.png) QQ-Plot when trying to fit a Normal distribution to my data ![alt text](http://i.stack.imgur.com/hUubx.png) **UPDATE:** From all the suggestions, I think I now understand why this cannot be a Cauchy distribution. Thanks to everyone. @HairyBeast suggested that I look at a negative binomial distribution so I plotted the following as well: QQ-Plot when Negative Binomial Distribution was used ![alt text](http://i.stack.imgur.com/DyHHV.png) Negative Binomial Distribution ![alt text](http://i.stack.imgur.com/DPNLZ.png)


title: Necessary and sufficient condition on joint MGF for independence

text: Suppose I have a joint moment generating function $M_{X,Y}(s,t)$ for a joint distribution with CDF $F_{X,Y}(x,y)$. Is $M_{X,Y}(s,t)=M_{X,Y}(s,0)⋅M_{X,Y}(0,t)$ both a necessary _and sufficient_ condition for independence of $X$ and $Y$? I seem to recall that it is, but I checked a couple of my textbooks and they only mentioned _necessity_ : $$F_{X,Y}(x,y)=F_X(x)\cdot F_Y(y) \implies M_{X,Y}(s,t)=M_X(s) \cdot M_Y(t)$$ That result is clear as independence implies $M_{X,Y}(s,t)=\mathbb{E}(e^{sX+tY})=\mathbb{E}(e^{sX}) \mathbb{E}(e^{tY})$. Since the MGFs of the marginals are determined by the joint MGF we have: $$X,Y\text{ independent} \implies M_{X,Y}(s,t)=M_{X,Y}(s,0)⋅M_{X,Y}(0,t)$$ But I found it hard to find a reference to the converse - or a counterexample! A bit of googling found an old Tennessee Tech page that states without proof or citation "it has been known that" the converse is true, which wasn't terribly reassuring! Is the following sketch proof workable? Given a joint MGF $M_{X,Y}(s,t)$, this uniquely determines the marginal distributions of $X$ and $Y$ and their MGFs, $M_X(s)=M_{X,Y}(s,0)$ and $M_Y(t)=M_{X,Y}(0,t)$. The marginals alone are compatible with many other possible joint distributions, and uniquely determine a joint distribution in which $X$ and $Y$ are independent, with CDF $F_{X,Y}^{\text{ind}}(x,y)=F_X(x) \cdot F_Y(y)$ and MGF: $$M_{X,Y}^{\text{ind}}(s,t) = M_X(s) \cdot M_Y(t) = M_{X,Y}(s,0)⋅M_{X,Y}(0,t)$$ So if we are given, for our original MGF, that $M_{X,Y}(s,t) = M_{X,Y}(s,0)⋅M_{X,Y}(0,t)$, this is sufficient to show $M_{X,Y}(s,t) = M_{X,Y}^{\text{ind}}(s,t)$. Then by the uniqeness of MGFs, our original joint distribution has $F_{X,Y}(x,y) = F_{X,Y}^{\text{ind}}(x,y) = F_X(x) \cdot F_Y(y)$ and $X$ and $Y$ are independent.


title: How do I check for bias of an estimator?

text: I need to check if an estimator $\hat\theta$ for the parameter $\theta$ is biased. Theory says I should compare the expected value of $\hat\theta$ versus the expected value of $\theta$. I assume the expected value of an estimator is the "weighted average" of the estimator and its distribution: $E[\hat\theta] = \int_0^\inf \hat\theta f(\hat\theta) d\hat\theta$*. If I'm right, to compute $E[\hat\theta]$ I need to know how $\hat\theta$ is distributed. For example: > $X$ is a random var with support $0 <= X <= \theta$ and pdf $f(x;\theta) = > 3x^2 / \theta^3$. Check if $\delta(x)=(4\bar X)/3$ is biased. *is it $E[\hat\theta] = \int_0^\inf \hat\theta f(\hat\theta) d\hat\theta$ or $E[\hat\theta] = \int_0^\inf \hat\theta f(x) dX$?


title: Confusion related to correlation in topic models

text: I was reading this paper related to correlated topic models. However, I didn't understand this correlation figure. In the figure given below, I am confused about the level curves shown in the simplex. In the first case, it says that there is diagonal covariance with non zero mean. In the second case, it says that there is negative correlation between components 1 and 2. In the third case it says that there is positive correlation between 1 and 2. I didn't understand how they intepreted it like this ![enter image description here](http://i.stack.imgur.com/yrx7P.png) In the figure given below, the topics are uncorrelated as it says ![enter image description here](http://i.stack.imgur.com/lzVlr.png) But I didn't get why it is so? Can anyone give me some insights?


title: How to calculate cumulative standards errors from several ones?

text: I want to estimate the cumulative effects of air pollution on mortality using GLM or Gam methods. With the following R code I was able to estimate the cumulative effect of PM10 on mortality which is estimated as 0.001513801 (1.51% increase in mortality for a 10 ug/m3 increase in PM10 ). My question is how can I build get the total standard error for estimates from lag 0 to 4 and build 95 % confidence interval for the cumulative estimate? Is it statistically valid to add up the standard errors just like the beta estimates from the output? I would also like to know if there are alternative methods of obtaining the total estimates with standard errors.               library(dlnm)      library(quantmod)     chic <- chicagoNMMAPS     fit <- glm(death ~ Lag(pm10, 0:4) + temp, data = chic,   family = poisson)     summ <- summary(fit)          Call:     glm(formula = death ~ Lag(pm10, 0:4) + temp, family = poisson,          data = chic)          Deviance Residuals:          Min       1Q   Median       3Q      Max       -4.0461  -0.8594  -0.0705   0.7390  22.4579            Coefficients:                           Estimate Std. Error  z value Pr(>|z|)         (Intercept)          4.749e+00  4.414e-03 1075.916  < 2e-16 ***     Lag(pm10, 0:4)Lag.0  8.395e-04  8.362e-05   10.040  < 2e-16 ***     Lag(pm10, 0:4)Lag.1  4.980e-04  8.993e-05    5.537 3.07e-08 ***     Lag(pm10, 0:4)Lag.2 -2.279e-05  9.007e-05   -0.253   0.8003         Lag(pm10, 0:4)Lag.3  2.198e-04  8.953e-05    2.455   0.0141 *       Lag(pm10, 0:4)Lag.4 -2.068e-05  8.312e-05   -0.249   0.8035         temp                -5.608e-03  1.520e-04  -36.902  < 2e-16 ***     ---     Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1          (Dispersion parameter for poisson family taken to be 1)      To obtain the total (cumulative ) effects for lags 0 to 4):               rown <- rownames(summ$coefficients)             tot <- grep("pm10", rown, fixed = TRUE)             coefs <- summ$coefficients[tot , "Estimate"]     total <- sum(coefs)     total      [1] 0.001513801      So in short how do I get total standard errors for lag 0:4?