title: Research design regarding similar names and comovements

text: I have a small question about a small research that I wanted to do regarding stock symbol similarities and their abnormal returns. I have about a sample of 10 000 stocks and I calculate how much they differ from each other (which gives a matrix of 10K x 10K with their "difference measure"). Now I was planning to test this by first grouping those stocks with a distance less than a certain threshold into a group and dropping those that have a larger threshold and then check whether these stocks have a similar return even though they are not in the same industry. I have a few questions, is this a good method to test for this? Or should I use the matrix of 10K x 10K? Because I need to check for totally unsimilar stocks too right? Furthermore am I forgetting something? I am trying to replicate the study done by Rashes (2001) on a larger scale.


title: Comparing density plots and scoring the combination

text: I have a set of density plots that contain the distribution of stock prices. Each graph has 5 density plots as follows that shows the distribution of monthly returns based on their ratings - a,b,c,d,e. A is the lowest rank while E is the highest rank. I need to come up with a score that tells me how different/away the 'A' plot is from the 'E' plot. i.e. ![Stock price distribution across different ranks](http://i.stack.imgur.com/IuwrX.png) Though the plots don't look very different, they are slightly different in the magnified version. I need to come up with a measure that can compare A's plot with E's plot. Although, mean, median, mode, skewness, kurtosis give the different aspects of a plot, I need one value that summarizes their differences. Any idea on what can be used? Thanks in advance!


title: Does it make sense to impute year of birth?

text: This is data cleaning and preparation stage question for me. I apologize if the question is basic, but I am a beginner. I have a dataset of a bit less than 4500 records. This is a survey and `year of birth` is an important field. Now 670 records do not have this information. I am inclined to think that I should treat this field as 'unknown' but I wanted to ask you a question: Does it ever make sense to impute year of birth? Perhaps you could also point me to any readings about whether demographic data can or should be imputed? Many thanks for your thoughts.


title: Difference between providing probabilities or providing a matrix to the chisq.test function in R

text: Take this example:               data <-matrix(c(227,751,193,541), ncol=2)     column1 <- c(227, 751)     probabilities <- c( 193/(193+541), 541/(193+541) )          chisq.test(data)     chisq.test(column1, p= probabilities)      when i apply the chi-squared test providing a matrix the results says that this is a > Pearson's Chi-squared test with Yates' continuity correction and provides a p-value of 0.158. when i perform the second chi-squared, providing the first column of the matrix and the probabilities calculated from the second column the both the results and the name of the test change dramatically: > Chi-squared test for given probabilities the reported p-value is 0.028. considering that i am trying to determine if the two datasets i have (columns in the matrix) are **NOT** different from each other: **what is the difference between these two tests? which one should i use?**


title: Replicator Neural Network, Step-wise function causing same prediction

text: In my project, one of my objectives is to find outliers in aeronautical engine data and chose to use the Replicator Neural Network to do so and read the following report on it (http://citeseerx.ist.psu.edu/viewdoc/download?doi=10.1.1.12.3366&rep=rep1&type=pdf) and am having a slight understanding issue with the step-wise function (page 4, figure 3) and the prediction values due to it. The explanation of a replicator neural network is best described in the above report but as a background the replicator neural network I have built works by having the same number of outputs as inputs and having 3 hidden layers with the following activation functions:   * Hidden layer 1 = tanh sigmoid S1(θ) = tanh,    * Hidden layer 2 = step-wise, S2(θ) = 1/2 + 1/(2(k − 1)) {summation each variable j} tanh[a3(θ −j/N)]    * Hidden Layer 3 = tanh sigmoid S1(θ) = tanh,    * Output Layer 4 = normal sigmoid S3(θ) = 1/1+e^-θ I have implemented the algorithm and it seems to be training (since the mean squared error decreases steadily during training). The only thing I don't understand is how the predictions are made when the middle layer with the step-wise activation function is applied since it causes the 3 middle nodes' activations to be become specific discrete values (e.g. my last activations on the 3 middle were 1.0, -1.0, 2.0 ) , this causes these values to be forward propagated and me getting very similar or exactly the same predictions every time. The section in the report on page 3-4 best describes the algorithm but i have no idea what i have to do to fix this, i don't have much time either :( Any help would be greatly appreciated. Thank you


title: Linear model in R - 2 fixed and one random factor

text: I have performed an experiment on plants and have some trouble analysing the data. For the experiment, we have three replicates (in time), two species and three treatments. In each experiment, we recorded the height of the plants. The dataset looks like this:                replicate | treatment | species | height      I assume `treatment` and `species` are fixed factors while `replicate` is a random one. For the analysis, I have used the procedure described here, so the code looks like this:               library(nlme)     rs <- read.table('mydata.txt', header=T)     rs$plot <- rs$species          m1.nlme = lme(height ~ species*traitement,               random = ~ 1|replicate,               data = rs)     summary(m1.nlme)      I got these results:               Linear mixed-effects model fit by REML      Data: rs          AIC      BIC   logLik       -326.2813 -310.415 169.1407          Random effects:      Formula: ~1 | replicate         (Intercept)   Residual     StdDev:  0.03603131 0.04329552          Fixed effects: height ~ species * traitement                                         Value  Std.Error  DF   t-value p-value     (Intercept)                       0.31546667 0.02225388 102 14.175806  0.0000     speciesmaize                      0.01645839 0.01127669 102  1.459506  0.1475     traitementtemoin                  0.21213868 0.01178502 102 18.000707  0.0000     speciesmaize:traitementtemoin    -0.04794414 0.01686908 102 -2.842131  0.0054      Correlation:                                   (Intr) gnF2F353 trtmnt     speciesmaize                     -0.249                     traitementtemoin                 -0.238  0.472              speciesmaize:traitementtemoin     0.166 -0.670   -0.702          Standardized Within-Group Residuals:        Min         Q1        Med         Q3        Max      -4.7398461 -0.3546214  0.1025950  0.6297446  1.6118499           Number of Observations: 108     Number of Groups: 3       The way I interpret this is that I have a `treatment` effect, as well as a `treatment*species` effect, but no effect of the `species` in itself. Is this correct?


title: Nonlinear regression: best transformation when getting very different parameter estimates

text: Disclaimer: Statistics is not my strong side, so if my question is nonsense I apologize. I'm a beginner, but really wanting to understand this. My question is: why do I get so widely different parameter estimates when using different transformations on my data in a non-linear regression ? I'm trying to do a nonlinear regression and to estimate the uncertainty of the fit (confidence interval) using linear approximation. From my understanding the more linear-like the shape of the nonlinear function, the more accurate will the confidence interval calculation by linear approximation be. I therefore want to transform the data to make it as linear as possible. The errors in $y$ can be assumed to be log-normal. My data is monotonic and assumed to follow a power function in most cases. $$ y = a*(x-x_0)^b $$ where $y$ is river discharge, $x$ is an arbitrary water level in the river and $x_0$ is the water level where where discharge $y$ is 0. This can be rewritten as log transformed, and nice and linear $$ log(y) = a + b \times log(x-x_0) $$. I need to estimate the parameters $a$, $b$ and $x_0$, so to do so simultaneously I use nonlinear regression. I also have some data that follows quadratic functions, so I would like to set up (and understand) a non-linear method. I use r and `nlsLM()` from `minpack.lm` to carry out the non-linear regression. Here is some example code:               library(minpack.lm)          xdata <- c(19,  21,  24,    25, 29, 34, 35, 40, 40, 46, 48, 48, 52, 56, 57, 65, 65, 68)     ydata <- c(10,  11, 14, 20, 24, 50, 42, 96, 89, 134,    135,    161,    171,    218,    261,    371,    347,    393)     df<-data.frame(x=xdata, y=ydata)          #weights applied in the case of no transformation (relative error assumed to be the same for all y data)     W<-1/ydata          # NLS regression with weights, no transformation     nlsmodel1<-nlsLM(y ~ a*(x-x0)^b,data=df,start=list(a=0.1, b=2.5,x0=0))          # log transformed     nlsmodel2<-nlsLM(log(y) ~ a+b*(log(x-x0)),data=df,start=list(a=0.1, b=2.5,x0=0))     > coef(nlsmodel1)               a           b          x0      0.005158377 2.719693093 4.896772931      > coef(nlsmodel2)             a         b        x0      -8.683758  3.445699 -4.139127           > exp(-8.683758)     [1] 0.0001693136      I understand that the weights are very important, and can have a say in the differences here, but not by this much? My judgement of the two parameter sets is that `nlsmodel1` performs "better", and that the `b` coefficient is too high in the fit from `nlsmodel2`. `nlsmodel2` does a poor job in the upper end of the data, with large residuals there. But why are they so different? I feel like I'm doing something very silly here, and is unable to see the error. I have tried some other transformations, for example only transforming LHS as `log(y)`, but the problem remains. I appreciate any tips that can help me improve, and not the least understand, the transformed fit. Cheers Related post #1 and post #2


title: Cardinal data as dependent variables

text: I wanted to find out if there are any implications to using OLS when i have cardinal data as dependent variables. So my dependent variables are counts of a certain outcome and they exist as natural numbers 0,1,2,3,...N My independent variable is a series of realized volatility. Is there a better way such data could be described?


title: Statistics for machine learning, papers to start?

text: I have a background in computer programming and elementary number theory, but no real statistics training, and have recently "discovered" that the amazing world of a whole range of techniques is actually a statistical world. It seems that matrix factorizations, matrix completion, high dimensional tensors, embeddings, density estimation, Bayesian inference, Markov partitions, eigenvector computation, PageRank are all highly statistical techniques, and that the machine learning algorithms that use such things, use a lot of statistics. My goal is to be able to read papers that discuss such things, and implement or create the algorithms, while understanding the notation, "proofs" and statistical arguments used. I guess the hardest thing is to follow all the proofs that involve matrices. What basic papers can get me started? Or a good textbook with exercises that are worth working through ? Specifically, some papers I would like to understand completely are :   1. Exact Matrix Completion via Convex Optimization, Candes, Recht, 2008   2. The Fast Cauchy Transform and Faster Robust Linear Regression, Clarkson et al, 2013   3. Random Projections for Support Vector Machines, Paul et al, 2013   4. High-Dimensional Probability Estimation with Deep Density Models, Rippel, Adams, 2013   5. Obtaining Error-Minimizing Estimates and Universal Entry-Wise Error Bounds for Low-Rank Matrix Completion, Király, Theran, 2013


title: Addition of multivariate gaussians

text: Given two multivariate gaussian (say in 2D with mean $\mu$ as a 2D point and convariance marix $\Sigma$ as $2$x$2$ Matrix) $N_1(\mu_1,\Sigma_1)$ and $N_2(\mu_2,\Sigma_1)$, I would like to derive the pdf of $N_1+N_2$. Can any one point me to the reference where i can find the pdf derivation of $N_1 + N_2$. Thanks in advance


title: Stratified sampling question

text: Good day all, Suppose that I am conducting a questionnaire study that is trying to measure level of awareness of subjects about a programming language and find the relation of those level of awareness to working conditions and methods etc. To improve my precision I decided to go with stratified sampling. If I have 1 criteria for stratification such as geo-distribution (to make sure I don't over represent subjects from areas that have less programmers), then I end up with 6 distinct strata (country provinces). I know how to go about analysing these to find margin error, standard error etc but I realised that is not good enough and I need to introduce more criteria's for stratification, such as level of education (so i don't over represent a group who are not very present between programmers based on their education), level of seniority etc. I have the proportion (%) for all these criteria's but I don't know how I go about sampling when I have more than one criteria? Thanks for your help :)


title: Additional features for regression

text: I want add to my regression additional features to complicate my model and lower bias, after searching around internet is seems a good idea is to add square features, that can help regression learn in coefficients more knowledge. I want add to my regression model square of pairs of features $x * y$, i think that can help learn my model relationships between features variables that can be highly correlated. Problem is when I see at interpretation of $x * y$ it seems incorrect, when i get $high * high = higher$ ex. $5 * 5 = 25$ that is good, but $low * low = higher$ seems strongly invalid ex. $-5 * -5 = 25$, should be $-25$ for my problem statement. Same problems I have for $high * low = low$, ex. $5 * -5 = -25$, i think should be somewhere between $<-5; 5>$, for example in $0$. That anyone have any thoughts about what equation I can use to get correct representation? I will appreciate other ideas for new virtual calculated features.


title: Binning raw data prior to building a logistic regression model

text: I have a set of variables for building credit scorecards with logistic- regression. I need to bin some variables, for e.g. years of credit history. What is the method to determine how many bins and what is the interval for each bin? Thanks


title: Statistical significance in time series and sub-series

text: So, when I did first year stats in undergrad, we did an experiment where we tampered with a bunch of coins, to see if it would cause a statistical difference in the results. This is a graph of the ratio of $heads:tosses$ for each series of flips: ![coin tosses](http://i.stack.imgur.com/0KGkp.png) We took the full data set in each case, and found that there were no significant differences from a null hypothesis of 50% heads. Had we stopped at 100 flips, or 150 flips we would have probably concluded that the Cupped coin was significantly biased. Would this have been invalid, and why? In particular, does it mean anything that the ratio is outside the 95% confidence interval more than 5% of the time?


title: Is a predictor with greater variance "better"?

text: I have a "basic statistics" concept question. As a student I would like to know if I'm thinking about this totally wrong and why, if so: Let's say I am hypothetically trying to look at the relationship between "anger management issues" and say divorce (yes/no) in a logistic regression and I have the option of using two different anger management scores -- both out of 100.   Score 1 comes from questionnaire rating instrument 1 and my other choice; score 2 comes from a different questionnaire. Hypothetically, we have reason to believe from previous work that anger management issues give rise to divorce.   If, in my sample of 500 people, the variance of score 1 is much higher than that of score 2, is there any reason to believe that score 1 would be a better score to use as a predictor of divorce based on its variance? To me, this instinctively seems right, but is it so?


title: What statistical test can I use to detect clumping?

text: I have time series data that represent dates/times of trades taken in a financial market. I would like to assign a score to this data that represents whether the trades are `mostly clustered` around particular time values or if they are `mostly spread out` evenly. I am going to have about 1000+ results per dataset. **Example situation one (High degree of "clustering" ):**               1. 01/01/01 : 13:00      2. 01/01/01 : 13:10     3. 01/01/01 : 13:15     4. 01/01/01 : 13:25     5. 03/05/01 : 17:20     6. 03/05/01 : 17:35     7. 03/05/01 : 17:40     8. 03/05/01 : 17:45      **Example situation two( Low degree of "clustering)"**               1. 01/01/01 : 13:00     2. 01/05/01 : 02:30     4. 02/12/01 : 06:40     5. 02/25/01 : 02:30     6. 03/30/01 : 21:10     7. 04/12/01 : 02:20     8. 05/02/01 : 03:25      I can of course convert all the timestamps to posix time or whatnot so doing calculation with the time values won't be a problem. I was thinking possibly standard error? (For those who want more background info: I am using backtest results to modulate the size of my entry position in a complex manner. If the results contain trades that are clustered together, then they don't really count as 1 trade each (more like one big trade). This means that such results are untrustworthy and I should not act on them.) Thanks!


title: Statistics of events and invitations

text: I am going to be hosting a number (~10) of potluck meals over the course of the summer, my pool of people to invite is about 40 people with about 10-15 coming to each meal. So I figure this would be a good opportunity to record data over time about the meals/people. The issue I am having is I am not sure what information to keep track of and what format to record it in. Here are some examples of trends I think would be interesting:   1. How many meals I have invited people to   2. On average which round of invites did people get invited to (some people rsvp as no in the beginning and so there is another 'round' of invites)    3. How many meals people have attended   4. What items people have brought I have started a spreadsheet where each page is a meal, the first few columns of the page represent different rounds of invites, I input a persons name in the column that corresponds to the round of their invite. The last two columns are the ultimate rsvp from any round of invitation and the item brought if applicable. To summarize I am looking for an efficient and concise way of recording the data associated with these meals for the trends mentioned. Additionally I am looking for other trends I can keep track of, I am doing a lot of this communication via email so timestamps would potentially be available for other interesting trends. Help with good tags for this question would be appreciated.


title: structural breaks in time series using matlab

text: in a plot of my time series there is clearly visible that there is structural break, but I have to find the exact date. I want test this with the chow test. Although I understand how to perform this test if the date of the structural break is know, by simply using a linear regression with two dummy's one for the intercept and one for the slope, $R_t$ = $\beta_0 $+ $\beta_0^* · D_i + \beta_1 R_{m,t} + \beta_1^* · D_t · R_{m,t} + ε_t,$ Than using the chow test.. But if I do not know the exact date, (in other words: I do not know when $D_i$ and $D_t$ should be 1) how can I find the exact date? Thank you very much for reply


title: Probability of x between two random variables

text: Given are $Z_1, Z_2$ i.i.d. standard normal. Find $P[Z_1 < t < Z_2]$ I have difficulties with working out how I should split the condition. Is $P[Z_1 < t < Z_2] = P[Z_1 < t, t < Z_2]$ or is some additional condition required? ## My question How should I approach this problem?


title: Calculate Hazard Ratio Manually for Continous Covariate

text: I found an example for manually calculating the hazard ratio when there are two groups.               d2 <- c(2, 2, 1, 2, 2, 0, 0, 4, 0, 2, 2, 0, 1, 0, 1, 1, 1)     d <- c(2, 2, 1, 2, 2, 3, 1, 4, 1, 2, 2, 1, 1, 1, 1, 2, 2)     r2 <- c(21, 19, 17, 16, 14, 12, 12, 12, 8, 8, 6, 4, 4, 3, 3, 2, 1)     r1 <- c(21, 21, 21, 21, 21, 21, 17, 16, 15, 13, 13, 12, 11, 11, 10, 7, 6)     l <- log(r2/r1)     d1 <- d - d2     dd <- cbind(d2, d1)     summary(glm(dd ~ offset(l), family = binomial))      where d is the number of deaths in total, and $d_i$ is the deaths in group $i$ and $r_j$ is the number of people in risk group $j$. I understand how this works. I would like to do a similar exercise for when there are not two groups of people, but only one vector of continuous data (e.g. expression of a particular gene) and no grouping factor. How can it be done for this scenario ?


title: Understanding an Implementation of normal conjugate in R

text: I am reading a code for a Bayesian clustering method. Both prior and the likelihood are normally distributed. If I understood correctly, such cases are called "conjugate priors" My question is about calculating the posterior mean and variance. So it is implemented in the code as following               d = nrow (data.i)##number of attributes     n = ncol (data.i)##number of replication for each attribute     Smatrix = matrix (1, nrow=d, ncol=d)##correlation of attributes     Imatrix = Smatrix - 1     diag (Imatrix) = 1          prior.precision = solve (sdWICluster^2 * Smatrix + sdTSampling^2 * Imatrix) #inverting the prior correlation matrix (prior.precision)     prior.mean = cluster.mean # mean of each clusters          sample.precision = sdResidual^(-2) * Imatrix     sample.mean = apply (data.i, 1, mean)#mean for each cluster          post.cov = solve (as.matrix(prior.precision + n*sample.precision)) # posterior covariance matrix     post.mean = as.vector (post.cov %*% (prior.precision %*% prior.mean + n*sample.precision %*% sample.mean)) # posterior of the mean      it seems the code has take the following formula, $\mu_{po} = C_{po}\times((\mu_{pr} \times \tau_{pr})+(n\times\tau_{li}\times\mu_{li}))$ $C_{po} = \tau_{pr} + n\times\tau_{li}$ As I said, this is the case of conjugates of normal distributions with unknown mean but variance known; however, it does not fit to the formula I have from here. (or may be it does by my eye is not able to catch). I appreciate if someone make some comments about the code and the formula


title: Plotting Options in R: Setting Axis Limits

text: I am plotting a density estimate for misclassification rate of some classifier using the standard `plot` and `lines` functions. Even though I've set `xlim=c(0.32,0.38)` and `ylim=c(0,100)` within the `plot` command, the x-limits and the y-limits are a bit wider than I'd like. How can I get the bottom left of the plot to have coordinates (0.32,0) and the top right of the plot to have coordinates (0.38, 100)?


title: Variance-covariance matrix for structural parameters in simultaneous equation models

text: The structural form of the linear simultaneous equations model simultaneous equations model can be written as $ \mathbf{y}_{i}^{\prime}\Gamma+\mathbf{x}_{i}^{\prime}\mathbf{B}=\epsilon_{i}^{\prime} $ which can be written in the reduced form model as $ \mathbf{y}_{i}^{\prime} =-\mathbf{x}_{i}^{\prime}\mathbf{B}\Gamma^{-1}+\epsilon_{i}^{\prime}\Gamma^{-1}\\ =\mathbf{x}_{i}^{\prime}\Pi+\mathbf{v}_{i}^{\prime}. $ The reduced form model can be estimated with different estimation methods and if the model is identified the structural coefficients $\widehat{\mathbf{B}}$ can be obtained by using the relation $\widehat{\mathbf{B}}=\widehat{\Pi}\widehat{\Gamma}$. I wonder how to drive the variance-covariance matrix of $\widehat{\mathbf{B}}$. I'd highly appreciate if you point out me any reference. Thanks in advance for your help and time.


title: Confidence interval for a proportion when sample proportion is almost 1 or 0

text: I know there are methods to calculate a confidence interval for a proportion to keep the limits within (0, 1), however a quick Google search lead me only to the standard calculation: $\hat{p} \pm 1.96*\sqrt\frac{\hat{p}(1-\hat{p})}{N}$. I also believe there is a way to calculate the exact confidence interval using the binomial distribution (example R code would be nice). I know I can use the prop.test function to get the interval but I'm interested in working through the calculation. Sample situations (N = number of trials, x = number of success):               N=40, x=40     N=40, x=39     N=20, x=0     N=20, x=1


title: Please help me with my MA(q) spreadsheet

text: Okay I have given this a go but I think it is totally wrong. Can someone please make this correct for me ? To understand the MA(q) model I really need to plug some numbers into it as I', finding these time series formula epically confusing as it is never explained what any of it actually is in the real world.                   Model : y = c + u + .5 et-1 + .25 et-2 + .125 et -3                   Ya  t   Yf t    mean U  e t          Actual Y at time t  Forecast Y at time t    1.69%   error at time t (et 0 is made up)             y = c + u + .5 et-1 + .25 et-2 + .125 et -3     Yf t - Mean u     0           1.69%   0.03     1   0.63%       1.69%   -1.06%     2   1.13%       1.69%   -0.56%     3   1.36%       1.69%   -0.33%     4   0.97%       1.69%   -0.72%     5   1.11%       1.69%   -0.58%     6   0.88%       1.69%   -0.81%     7   0.90%       1.69%   -0.79%     8   1.06%       1.69%   -0.63%     9   1.38%       1.69%   -0.31%     10  0.96%       1.69%   -0.73%     11  1.29%       1.69%   -0.40%     12  1.39%       1.69%   -0.30%     13  0.98%       1.69%   -0.71%     14  1.42%       1.69%   -0.27%     15  1.21%       1.69%   -0.48%     16  1.04%       1.69%   -0.65%     17  1.45%       1.69%   -0.24%     18  1.11%       1.69%   -0.58%     19  1.92%       1.69%   0.23%     20  1.89%       1.69%   0.20%     21  2.00%       1.69%   0.31%     22  1.40%       1.69%   -0.29%     23  2.15%       1.69%   0.46%     24  2.11%       1.69%   0.42%     25  1.76%       1.69%   0.07%     26  2.55%       1.69%   0.86%     27  1.80%       1.69%   0.11%     28  2.23%       1.69%   0.54%     29  2.39%       1.69%   0.71%     30  1.55%       1.69%   -0.14%     31  2.67%       1.69%   0.99%     32  2.50%       1.69%   0.81%     33  1.88%       1.69%   0.20%     34  2.18%       1.69%   0.49%     35  2.42%       1.69%   0.73%     36  1.75%       1.69%   0.06%     37  2.54%       1.69%   0.85%     38  2.03%       1.69%   0.34%     39  2.08%       1.69%   0.39%     40  2.58%       1.69%   0.89%     41  1.71%       1.69%   0.02%     42  1.56%       1.69%   -0.13%     43  1.98%       1.69%   0.29%     44  1.85%       1.69%   0.16%     45  1.66%       1.69%   -0.03%     46  2.19%       1.69%   0.50%     47  1.52%       1.69%   -0.17%     48  1.75%       1.69%   0.06%     49  1.92%       1.69%   0.23%      FORECAST 50 1.80% 1.69% 0.11% FORECAST 51 1.81% 1.69% 0.12% FORECAST 52 1.81% 1.69% 0.12% FORECAST 53 1.79% 1.69% 0.10% FORECAST 54 1.79% 1.69% 0.10% FORECAST 55 1.78% 1.69% 0.09% FORECAST 56 1.77% 1.69% 0.08% FORECAST 57 1.76% 1.69% 0.07% FORECAST 58 1.76% 1.69% 0.07% FORECAST 59 1.75% 1.69% 0.06% FORECAST 60 1.75% 1.69% 0.06% FORECAST 61 1.74% 1.69% 0.05% I have tried to do this to the step that I am upto before estimating the first error term for t = 0 and parameters of the model but it does not look right. Hope you can follow this it is hard to get the formatting right. It just looks like it ends up being too similar to the mean to be useful. I want something which tries to predict the random patterns seen in the time series, for example you can see when the point goes up, it is very rarely higher in the very next point. The formula for the forecasting is basically to take the mean and multiply it by the last 3 error terms and multiply by the coefficients that I just made up. Thanks. Thanks. As a side note I really find the notations and explanation about this to be very confusing. Thanks.


title: Rejection sampling from a normal distribution

text: I am running a Monte-Carlo simulation and I sample from various normal distributions. I was just wondering, is there a way by which I can increase the probability of selecting a point from the tails (i.e., $[-5\sigma,-3\sigma]$ and $[3\sigma,5\sigma]$) as compared to the probability associated with the interval $(-3\sigma, 3\sigma)$? Edit: ![enter image description here](http://i.stack.imgur.com/UzAEY.jpg)


title: How do I handle predictor variables from different distributions in logistic regression?

text: I am using logistic regression to predict y given x1 and x2:               z = B0 + B1 * x1 + B2 * x2     y = e^z / (e^z + 1)      How is logistic regression supposed to handle cases in which my variables have very different scales? Do people ever build logistic regression models with higher-order coefficients for variables? I'm imagining something like this (for two variables):               z = B0 + B1 * x1 + B2 * x1^2 + B3 * x2 + B4 * x2^2      Alternatively, is the right answer to simply normalize, standardize or rescale the x1 and x2 values before using logistic regression?


title: Why is it not OK to do a Pearson correlation on proportion data?

text: An online module I am studying states that one should **never** use Pearson correlation with proportion data. Why not? Or, if it is sometimes OK or always OK, why?


title: Similarity between two samples of spacialized data ignoring local fluctuations

text: I have data from a set of sensors, and would like to measure similarity between two samples of data ignoring local fluctuations. My problem is, my sensors are not regularly placed : some are quite alone, others are many in the same area, there is no regular pattern in their locations. I want to measure similarity in terms of large scale, local fluctuations (even of high amplitude) are much less important as long as they are balanced by other nearby fluctuations. I though about computing average values of a point and its neighbours, and making my indicator as the sum of absolute value of differences between the local averages of my two samples, but I have not idea whether it's a good approach or not.


title: Is p-value also the false discovery rate?

text: In http://surveyanalysis.org/wiki/Multiple_Comparisons_(Post_Hoc_Testing) It states "For example, if we have a p-value of 0.05 and we conclude it is significant the probability of a false discovery is, by definition, 0.05." My question: I always thought false discovery is Type I error, which is equal to the chosen significance levels in most tests. p-value is the value calculated from the samples. Indeed Wikipedia states "The p-value should not be confused with the significance level α in the Neyman–Pearson approach or the Type I error rate [false positive rate]" So why the statistician in the link claims Type I error is also the p-value?


title: Correlation versus cause-effect regression

text: I know correlation does not imply causation. I have read it nth time. (i.e. weight does not cause height etc. etc.) However, to find the effect of a moderator variable on X-Y relationship, a regression model is used such as the GLM in SPSS to test for interaction or multiple regression. See here. **My Question: If the X-Y relationship is a correlation, then why is a cause- effect model used in this instance?** As far as I understand, it makes little sense to classify a variable as independent or dependent in a correlation analysis. I apologise if this seems like a 'silly' question. In my past life, I had often told my students that there are no silly questions; just questions!


title: Why is it futile to use the deviance as a goodness-of-fit measure for Bernoulli data?

text: In Ordinal Data Modelling by Johson & Albert, page 102-103: > For Bernoulli observations [...] the asymptotic chi-squared distribution of > the deviance statistic may not pertain. Indeed, for the linear logistic > models with Bernoulli observations, the deviance function can be expresses > solely as a function of the MLE of the regression parameter, **which > demonstrates the futility of using this statistic to measure goodness of fit > for such data**. Could someone kindly explain what this means, please?


title: Estimation with transformations of variables after multiple imputation

text: I would really appreciate some advice. I have two sequential measurements of severity of illness $(s_1, s_2)$ that are incompletely observed along with an outcome measure ($y$). I am trying to estimate if the change in the measurement $(d = s_2 - s_1)$ adds anything to a model based on the second measurement $s_2$ alone. There is an extra complication in that the severity of illness measurement is a score derived by categorising 10 continuous measurements of physiology ($x_1,x_2...x_{10}$), and weighting them. This is helpful because it is the standard clinical method in which to report severity, and it handles variables that are not always linearly associated with survival (e.g. very high or very low blood pressure is a bad thing --- you're best off in the middle). I have been trying to use multiple imputation assuming MAR in Stata 12. The final model would therefore be $y=\alpha + \beta_1s_2+\beta_2d+e$. However, I think there are several different imputation strategies and I am not sure which is correct.   1. I use the raw physiology in the imputation model, and then derive the weights afterwards.   2. I derive the individual weights ($x_1,x_2...x_{10}$) first, and then use these in the imputation model.    3. I use the aggregate score in the imputation model.  I am aware that using a transformation of a variable in the final model after the imputation (e.g. $x^2$ when you only included $x$ in the imputation) biases the effect of $x^2$ to zero. I am not sure if this is also the case with a categorisation and summation? In which case, I assume (2) would be the best option. Otherwise (1) seems best to me. I have already discounted (3) because it seems to lose a great deal of the richness of the data. Thanks in advance.


title: Ranks test problem

text: This is a questions I got for homework in _test planning and variance analysis_. I apologize in advance for the way my professor phrases things :D I would love to hear hints or tips on how to approach this. So here goes: > A statistician decides to perform the next a-parametric test to examine the > $H_0$ of variance Analysis. First he will calculate the test statistic in > the following way: Calculating the ranks of all observations with no > affiliation to a group, and calculating F statistics of the anova test on > the ranks instead of on the observations. Then, he will calculate the > p-value using a permutation test, since under $H_0$ there is the same chance > for all the permutations of ranks when dividing into groups. > >   1. What is the SST calculated on the ranks of observations? Show that the > SST does not depend on the observations, only on the number of observations. > >   2. Prove that the p-value of this test will be equal to the p-value that > is calculated with permutation on the Kruskal-Wallis test statistic. Clue: > Show that the test statistic of Kuskal-wallis is SSB/(SST/N-1) when SSB and > SST are calculated on the ranks of observation and N is the total number of > observations. > > This is what I've got for (1). Not sure about it but it makes sense to me.. What I don't get is how to show that SST does not depend on the observations, because obviously we need them in order to calculate the SST. $k$ - number of groups   $n_i$ - number of observations in group $i$   $\bar{d}$ - mean of ranks of all groups   $\bar{d}_i$ - mean of ranks in group $i$ $$\operatorname{SSB} = \sum_{i=1}^k n_i(\bar{d}_i-\bar{d})^2$$ $$\operatorname{SSW} = \sum_{i=1}^k \sum_{j=1}^{n_i} (d_{ij}-\bar{d}_i)^2$$ $$\operatorname{SST} = \sum_{i=1}^k \sum_{j=1}^{n_i} (d_{ij}-\bar{d})^2$$ And actually when I think about it, I'm not so sure what is difference between the test that is described in this question and Kruskal Wallis.


title: Discrete Fourier Transform and uneven sampling

text: In this blog article an example is given of one can use the DFT to detect frequencies much higher than the sample rate. In the comments sections I asked how it was done, since DFT normally requires evenly sampled data. The author of the blog responded concisely with: > One can still define the DFT as proportional to > > $$ X(\nu) = {1 \over N} \sum x_n e^{-i 2 \pi \nu t_n} $$ How would one do that? I am asking for advice to do exactly that. I can replicate the example if I use the data compensated method of Ferraz-Mello 1 but the blog author says that normally DFT will work using that hint. 1 http://articles.adsabs.harvard.edu/cgi-bin/nph- iarticle_query?1981AJ.....86..619F&data_type=PDF_HIGH&whole_paper=YES&type=PRINTER&filetype=.pdf


title: Explanation of minimum observations for multiple regression

text: I feel like every question I've asked on CrossValidated has lead back to looking at the number of observations I have per variable. I understand that there are many rules of thumb out there depending on your field, your expected effect size, etc. These come back to recognizing that too few observations can lead to overfitting. I understand how overfitting can be a problem when looking at a regression model with a single predictor. Having only two observations leads to a perfect answer, whereas using least squares to solve an overdetermined model leads to generalization. However, I have trouble making sense of how overfitting can still be a problem if you have 20 predictors and, say, 31 observations. It seems like you have mitigated the problem with 10 extra observations, but I suspect I am missing something in how least squares solves overdetermined systems. What I assume follows if it is true that a system is overfitted is that the relationships between predictors explained by the betas also do not hold. Finally, if overfitting is a problem by having not enough observations, can this be solved by using forward stepwise regression? Or is there a good possibility of missing significant predictors due to eventually reaching a point where there are too many predictors and not enough observations? To recap:   1. Can someone explain why the case where the number of predictors is $N-1$ where $N$ is decently large (say, $N > 10$) can still overfit? Why is it not that you need a minimum number of observations **regardless** of the number of predictors?   2. Can overfitting be solved by using forward stepwise regression?


title: Simple Logistic Regression Questions

text: # Problem I want to estimate the probability of a user choosing a car or transit as a transport mode. Following is the output of the simple logistic regression: ![enter image description here](http://i.stack.imgur.com/9qQZC.jpg) I want to estimate the probability of choosing car (auto) if INVT for auto = 3, INVT for transit = 7, OVT for auto = 1, OVT for transit = 2.5, OPC for auto =2 and transit fare = 3. # What I tried For the probability of choosing car when the user owns a car and work in suburbs, I am using INVT for auto, OVT for auto, OPC for auto and transit fare along with AO as 1 and DW as 0. > The probability = Pr = exp(1.45 - (0.00897*3) -(0.0308*1) - (0.0115*2) > -(0.00708*3) +(0.77 * 1) - (0.561*0))/(1+exp(1.45 - (0.00897*3) -(0.0308*1) > - (0.0115*2) -(0.00708*3) +(0.77 * 1) - (0.561*0))) > > Pr = 0.893 # Question Is this the right approach to find this probability? Or should I also include the INVT, OVT and OPC for transit as well?


title: Design a feature with time and presence information

text: Context: I am working on a decision tree classifier, trying to classify businesses as to whether they are likely to have an event occur (default) in the next 90 days. One input I get is whether, and when, the business is scheduled to have a related event (let's call it "warning"). A business may have no warning scheduled, or it may have a warning scheduled some number of days in the future (e.g., we learn today that a warning has been scheduled for 15 days hence). If a warning is scheduled further out, we may presume that the business has more time to cope and fix it, so it is "worse" to have 15 days till warning, than it is to have 85 days till warning. I'm looking for a way I can massage/transform this data into a single, real- valued feature for use in a decision tree, that 1. preserves the natural ordering that intuitively exists here, and 2. won't make for too ugly of statistics? Possibilities:   1. Just use the number of days until scheduled warning as a feature. Pros: preserves order (except for no warning) and natural interpretability. Cons: when "no warning" exists, "zero" is a particularly poor value to use, since it violates the ordering.   2. Use "+infinity" or some arbitrarily large number of days when "no warning" exists, but otherwise use number of days. Pros: preserves order and interpretability (except for no warning). Cons: most values are going to be in the 1-90 range, so choosing a representative value for "no warning" skews the data.   3. Use "zero" for no warning, and transform the number of days till warning with the inverse. Pros: smooth, well-ordered over the domain, output range of 0,1, "no warning" is interpreted as "zero" effect. Cons: loses native interpretability of values, unsuitable for other (linear?) models.   4. Use two separate features, one a binary for warning, one the days till warning. Pros: properly separates the concerns (?). Cons: still doesn't deal with the question of what to put in for missing value in days till warning. If, as I suspect, this is something of a rookie problem, any pointers for standard approaches to this kind of problem are welcome, including "RTFM" responses, though specific sources, technique names, or search terms would be quite welcome.


title: Statistical thresholds for large time-series

text: I am working with matrices, of the size N*M, where each cell corresponds to the Pearson's correlation between two time series. I want to threshold each matrix such that it would retain only significant correlations. Now, my problem is that the time series are very large (>1500 data points) and thus rather weak correlations turn out to be (very) significant. I am already considering only correlations whose p val survives correction for multiple comparison (to correct for the number of tests done within each matrix), however this still leaves quite a few weak correlations in each matrix. I prefer not to choose a threshold arbitrarily, so I was wondering if anyone could direct me to established methods for determining statistical thresholds when working with very large time series Many thanks for your help!


title: Can neural net extrapolate output value

text: Let's assume we have a training set with $y \in \mathbb{R}$. Thus all the data is between $y_{min}$ and $y_{max}$. If we built a decision tree model it cannot return $y_{pred}$ outside the given range (using any combination of input features). Thus decision tree cannot extrapolate in terms of predicted values. Can a neural net regression model extrapolate and return $y_{pred}$ values outside the $y$ range in a training set? Does it depend on the activation function or not? Below is my attempt to answer this question. The output neuron of the model is just $\sum \Theta_ia_i$, where $\Theta_i$ - weight of i-th neuron on the previous hidden layer, and $a_i$ - value of activation function of that neuron. If we use logistic function then $a \in (-1;1)$. Thus maximum possible $y_{pred} = \sum \Theta_i$, assuming that all $a$ reach their maximum value around 1. But if we will use linear activation function, which doesn't have restrictions on output values of $a$ ($a \in \mathbb{R}$) the model will return $y_{pred} \in \mathbb{R}$, which can be ouside $y$ range of the training set. Is my line of reasoning correct or there are some mistakes?


title: Endogenous variable and statistical significance in OLS?

text: I was wondering what the following OLS scenario would imply: a variable is endogenous (i.e. correlated with the error term) yet is statistically significant. Alternatively, what if in, once again OLS, the variable is exogenous yet statistically insignificant (although the overall regression is indeed statistically significant - by say an F test). Thanks!


title: How do i interpret Johansen's test? (Matlab)

text: I found several similar questions with no answer. But I'm gonna give it a shot anyway. I'm trying to find assets (stocks) which are co-integrated, since Engle- Granger's test may give cause to spurious regression it's suggested by several articles that you use Johansen's test as well to ensure a long-term co- integration. Even though I've not taken any course in time series analysis yet the Engle-Granger's test was not to difficult to understand. Johansen's test was on the other hand slightly more challenging, so I'm trying to understand how the interpret the results you're given in matlab. So when I run the test for a nx2 vector in Matlab I get a 1x2 matrix as result with the inputs true or false. So the question is when I can reject/not reject the hypothesis.               ************************      Results Summary (Test 1) Data: Effective sample size: 13 Model: H1 Lags: 0 Statistic: trace Significance level: 0.05 # r h stat cValue pValue eigVal 0 1 23.6646 15.4948 0.0031 0.7486   1 1 5.7169 3.8415 0.0170 0.3558 h =                     r0       r1        t1    true     true       pValue =                     r0           r1           t1    0.0031234    0.016977      stat =                     r0        r1         t1    23.665    5.7169      cValue =                     r0        r1         t1    15.495    3.8415      mles =                     r0              r1               t1    [1x1 struct]    [1x1 struct]


title: MPNet output to visualization

text: I've been using igraph and statnet to do most of my network analysis, but I am trying to replicate results that I use MPNet (an expansion to PNet) in order to perform some simulations on a multi-level network. I'm hoping someone has experience with PNet or MPNet and can offer a suggestion. I've run a few simulations of small networks, and the output for any specific simulation doesn't really lend itself to being visualized easily in igraph and statnet outside of manually inputting the networks and visual attributes of each vertex. This is problematic for me because I want to simulate large-ish (n>20 in each level) networks. _However,_ it does appear that the output is intended to be visualized by the way it is written (see below), but I'm unsure what program it's intended for and cannot find much existing MPNet (or PNet) help online. **My question is this:** is there a way to use MPNet output for a specific simulation to visualize the network that is being simulated without manually entering in each network and visual attributes into some package? Or in other words, did MPNet designers have a destination in mind with their output? * * * The following picture is a quick MSPaint rendition of this 10 vertex multi- level simulated network (because we all need a little MSPaint in our lives): ![enter image description here](http://i.stack.imgur.com/Lts6N.jpg) I have four output text files for a specific simulation (sorry in advance for the huge amount of output): Network A:                 *vertices 5     1 "" box ic Blue bc Black     2 "" box ic Blue bc Black     3 "" box ic Blue bc Black     4 "" box ic Blue bc Black     5 "" box ic Blue bc Black     *matrix     0 1 1 1 1      1 0 0 0 0      1 0 0 1 1      1 0 1 0 1      1 0 1 1 0           Density 0.7000     Degree Distribution(s)     Mean degree 1.4000     Stddev  1.9105     Skewness    0.8547     Clustering coefficent(s)     Global  0.8000      Network B:               *vertices 5     1 "" ellipse ic Red bc Black     2 "" ellipse ic Red bc Black     3 "" ellipse ic Red bc Black     4 "" ellipse ic Red bc Black     5 "" ellipse ic Red bc Black     *matrix     0 0 1 0 0      0 0 1 0 1      1 1 0 1 1      0 0 1 0 0      0 1 1 0 0           Density 0.5000     Degree Distribution(s)     Mean degree 1.0000     Stddev  1.6583     Skewness    1.2718     Clustering coefficent(s)     Global  0.3750      Network X:               *vertices 10 5     1 "" box ic Blue bc Black     2 "" box ic Blue bc Black     3 "" box ic Blue bc Black     4 "" box ic Blue bc Black     5 "" box ic Blue bc Black     6 "" ellipse ic Red bc Black     7 "" ellipse ic Red bc Black     8 "" ellipse ic Red bc Black     9 "" ellipse ic Red bc Black     10 "" ellipse ic Red bc Black     *edges     1 6     1 9     3 7     3 8     3 9     3 10     4 6     4 7     4 8     5 8     5 9          Density 0.4400     Degree Distribution(s)     Mean degree A   2.2000     Stddev A    1.4832     Skewness A  -0.2648     Mean degree B   2.2000     Stddev B    0.8367     Skewness B  -0.2459     Clustering coefficent     Global  0.2667      Network M:               *vertices 10     1 "" box ic Blue bc Black     2 "" box ic Blue bc Black     3 "" box ic Blue bc Black     4 "" box ic Blue bc Black     5 "" box ic Blue bc Black     6 "" ellipse ic Red bc Black     7 "" ellipse ic Red bc Black     8 "" ellipse ic Red bc Black     9 "" ellipse ic Red bc Black     10 "" ellipse ic Red bc Black     *edges     1 2 1 c Blue     1 3 1 c Blue     1 4 1 c Blue     1 5 1 c Blue     3 4 1 c Blue     3 5 1 c Blue     4 5 1 c Blue     6 8 1 c Red     7 8 1 c Red     7 10 1 c Red     8 9 1 c Red     8 10 1 c Red     1 6 1 c Black     1 9 1 c Black     3 7 1 c Black     3 8 1 c Black     3 9 1 c Black     3 10 1 c Black     4 6 1 c Black     4 7 1 c Black     4 8 1 c Black     5 8 1 c Black     5 9 1 c Black      *Also, if someone has enough reputation and feels the need to add the tags PNet and/or MPNet to this post so that it could be more specific, I would appreciate it.


title: Statistical Difference from Zero

text: I have a set of data that represents periodic readings. The data shows an upward trend but I need to test for a statistical difference from zero. I believe I should use the t test, two tailed but what should I use for the second set of data? Zeros, the starting value?               0.2245     0.243     0.2312     0.1795     0.1923     0.17     0.2025     0.2059     0.2394     0.205     0.2201     0.2261     0.1817     0.2143     0.2126     0.237     0.1984     0.228     0.2292     0.2236     0.2096     0.2258     0.2155


title: better precision worse recall

text: is it possible for an algorithm A to have a better precision but worse recall (or better recall but worse precision) than another algorithm B? Although I know that precision and recall are different things, it seems to me that if algorithm A has a better precision (or recall) than algorithm B then algorithm A will also have a better recall (or precision) than B. Thanks Ahmet


title: Probability of observing a measurement with error under a non-normal distribution

text: Suppose I have some kind of distribution which gives the probability of observing a given value (assume that this distribution is **not** normal). Now, I have a certain observation and I want to calculate how likely it is to observe a value that is equally or more extreme under this distribution. So, I am thinking of doing this empirically by turning my observation into a quantile and calculating the probability empirically under the distribution of: quantile_distribution < quantile_observation and quantile_distribution > abs(1 - quantile_observation) I'm not sure whether that is a valid way to do it, but for the second part I'm going to assume it is. Now, here's where it gets a bit more complicated: My observation has some sort of measurement error associated with it. I can model this error with a probability distribution (this time it is roughly normal). Given this spread, how would I calculate a p-value as above (or would it be a distribution of p-values, which I have never seen)? I've been hacking something together where I take the p-value of observing the mean and the 95% CI values under the first distribution and then just take the maximum and report that, but obviously there has to be a better way to do it


How can I find the standard deviation in categorical distribution, where the elements have non-numerical attributes (e.g. colors)? For example, I have a bag of marbles with $n$ colors. There's an infinite number of marbles in the bag, and are biased towards a certain color with a probability of $\frac{x}n, x>1$. From the bag, I pick $m$ marbles and get the probability distribution with respect to their colors. From this distribution, by picking $q$ colors with the highest probabilities, I want to convince others with a confidence interval $k$% (e.g. 95%), that one of the colors I chose is the one that the marbles are biased towards. In this scenario, what are some of the analysis techniques that can be used to find $q$, given $x,n,k,$ and $m$?


title: Which statistics textbooks would explain why linear least squares regression fails when data is heteroskedastic

text: I have recently interviewed for a statistical analysis job and was asked a question about why linear least squares regression fails when the data is heteroskedastic. The correct answer to this question, according to the interviewers, is that heteroskedastic data means that the equation of the regression line produced by least squares regression is an unbiased estimator of the true relationship, but that it is NOT efficient, essentially because the part of the dataset where the variance is smaller than average is effectively underweighted. My question is **which textbooks could I use to find more detail about this topic, and other similar topics at this level** , e.g.   * the relationship between data being normally distributed and least-squares linear regression being the maximum likelihood estimator for the straight line fit  [I have a degree in mathematics but with minimal statistics background & understand general probability concepts such as the central limit theorem, random variables, etc, and I know high school level statistics up to the British A-level S4 statistics, however I lack a certain level of statistics knowledge and don't know what I don't know or where to find out more... ]


title: Help in convolution operation

text: I am stuck in understanding how the output of convolution operation is obtained. Can somebody please show the steps? The question is an image is processed by applying 3*3 mean filter . What is the impulse response? The formula is $g(m,n) = \sum_{\alpha = 0}^{M-1} \sum_{\beta=0}^{N-1} f(\alpha,\beta)h(m-\alpha,n-\beta)$ The image function $f$ =                                        0 0 0 0 0                              0 1 1 1 0                              0 1 1 1 0                              0 1 1 1 0                              0 0 0 0 0          h = 1 1 1         1 1 1         1 1 1      The result of the operation (obtained in Matlab) is g(m,n) =               0 0 0 0 0 0 0     0 1 2 3 2 1 0     0 2 4 6 4 2 0     0 3 6 9 6 3 0     0 2 4 6 4 2 0     0 1 2 3 2 1 0     0 0 0 0 0 0 0      I cannot understand the steps how to do it in paper, can somebody please show the few steps how to go about it ? Based on the reply and the formula, I have worked out one example, but the answers don't match with the formula and the other technique. For convolution at pixel f(0,0) the answer = 0. But I am getting 0 using the formula and 1 using the drag method. The cell to which this result will be mapped is g(0,0) and I have put the answer there. For convolution of pixel element at f(1,1) = 1 also, when I evaluated I am getting answer 2 instead of 4. The result should go to g(2,2). But how do I get values for g(5,0) since the original image f has no fifth row!! How to calculate values for fifth and sixth rows? Please let me know what is wrong in my understanding. Is my coordinate system incorrect? Shall really appreciate. ![enter image description here](http://i.stack.imgur.com/FKado.jpg)


title: "Explaining" a time series - conceptual explanation needed

text: Imagine I have a time series for an animal population and a time series for a climatic variable during the same time period and at the same location. Unfortunately the data are observational (i.e., no experimental manipulation). Now imagine there is a decreasing linear trend in both variables (and no seasonality). I know before hand that these variables will be correlated because they share a time trend (e.g.,http://stats.stackexchange.com/a/8037/38125). But is there anything more I can infer from this data? I am looking for a rather low-level conceptual explanation as to whether this is possible and why. I have 2 thoughts on how this might be accomplished, but haven't figured on my own if either would be valid.... 1) Regress population against both time and the climate variable, or 2) "Detrend" both 'population' and 'climate variable' and regress their residuals against each other I can provide more detail if necessary but wanted to keep this original post brief. Thanks for your help!


title: R Interpret coefficient from Survreg(dist="gaussian")

text: I was wondering if anyone could help me on how to interpret the coefficient from an analysis I have carried out in R (`survival` package). The data is right censored, the dependent variable (time to event) was modelled under the gaussian distribution and the independent variable is a categorical variable (3 categories). my output is as follows:               coeff   se     p-value     -0.107  0.048  2.7x10-2      Can anyone tell me how to interpret this coefficient under the following model?


title: nlme regression & curve comparison in R. ANOVA p-value

text: I obtained data of typical time point in two conditions in this shape:               weight   condition   time     0.1307857    Transf    1     0.1926429    Transf    2     0.2734286    Transf    3     0.4403571    Transf    4     0.6037143    Transf    5     0.9036429    Transf    6     1.5454286    Transf    7     0.1370714       Unt    1     0.2005000       Unt    2     0.2973571       Unt    3     0.4592143       Unt    4     0.8336429       Unt    5     1.3099286       Unt    6     2.1470000       Unt    7      I am using the package `nlme` in `R` to test the difference between the two conditions in time               fm1 <- lme(weigth~condition,random=~1|time,data=data_new)     anova.lme(fm1, adjustSigma = F)      I got this:               numDF denDF  F-value p-value     (Intercept)     1     6 8.421439  0.0273     condition       1     6 4.208786  0.0861      Questions:   1. Am I doing things correctly with this procedure to test the differences between the two conditions in time?   2. If yes, which is my p-value, the intercept one or the condition?   3. If this is not the proper test, what are the alternatives?


title: Comparison of differences between pairs of samples of unequal size

text: I have four samples $x_1, x_2$ and $y_1, y_2$ with $n_{x1} \neq n_{x2} \neq n_{y1} \neq n_{y2} $. I calculated, using a Wilcoxon rank sum test, that $x_1$ is significantly different to $x_2$ and $y_1$ significantly different to $y_2$. However, I would like to test whether the difference $x_1 - x_2$ differs significantly to $y_1 - y_2$ but I have no idea how do to that given the unequal sample sizes. Any ideas or suggestions would be really appreciated.


title: A stopping rule for permutation test

text: I perform a permutation test multiple times on different datasets, each time I am only concerned about significant $p$ values. To reduce computation time would it be correct to introduce this kind of stopping rule: After a certain number of $N$ permutations to check whether $p$ is greater than a particular value. So, for example if $p>0.1$ after $N$=200 permutations then the lower bound of 95% confidence interval would be greater than $0.05$. Therefore, calculations could be stopped as a true $p$ is not significant. Just want to make sure I am doing it right. Thank you.


title: Finding latent variables in 7 dichotomous variables

text: I have read other posts on conducting Factor analysis (FA) with dichotomous variables and although it appears clear that FA done in the default way is not appropriate, I am still unclear about a few things. I have 7 variables that are yes/no responses. They are assessing peoples functioning (e.g., working yes/no, involvement with justice system yes/no etc). I want to combine these responses to have two variables that represent domains of functioning (e.g., occupational and social). I suspect three of them will be occupational and 4 of them social. I would like to do something analogous to a FA to confirm this. I only have access to SPSS (and am only familiar with this package) so am looking for advice on how to do this with SPSS. Do I do a FA with tetrachoric correlations instead of Pearsons (how is this done?)? Or can I do Latent Trait Analysis in SPSS? I am specifically looking for references with step-by-step instructions.


title: How to form a confidence band around the trend fitted from time series data

text: I have a time series data set. I can decompose it and get the trend but I would like to put confidence ranges around the trend (past) not the forecast- ed component. The decompose function also doesn't handle N/As very well so is there another way to define the trend with data that has N/As. I have been trying to use Holt-Winters and ARIMA but neither seem to be able to do this.


title: Fitting the parameters of a stable distribution

text: I have a data set and I have to fit this data set with a stable distribution. The problem is that the stable distributions are known analytically only in the form of the characteristic function (Fourier transform). How can I do this?


title: How to extract daily values from 16-day dataset throughout a year?

text: I have two sets of data. The one is daily temperature data and another is 16-day vegetation observation data (value of Jan 1, Jan 17 etc), throughout a year. But I need daily values of vegetation observation data. How to obtain daily value from 16-day vegetation data?


title: How to find which variables are most correlated with the first principal component?

text: I came across an article where the authors did a Principal Component Analysis on gene expression data, and found out the genes that are most correlated to the 1st principal component, and they used that gene list for further analysis. Can somebody tell me how to find out entities (genes in this case) that are most correlated to the 1st principal component? Here's the link to the original free article, and this is how they've calculated it: > Results of the gene set enrichment analysis of the genes most correlated to > the 1st principal component. The correlations of the genes with the 1st > principal component were transformed to SDs from the mean, then genes with > values > 1.5 (positive correlation) or < -1.5 (negative correlation) were > selected. I can make a toy example of 6 samples and 20 gene matrix and do the PCA the following way but how to proceed next:               rm(list=ls())     set.seed(12345)     my.mat <- matrix(rnorm(120,0,0.5),nrow=6,byrow=TRUE)     rownames(my.mat) <- paste("s",1:6,sep="")     colnames(my.mat) <- paste("g",1:20,sep="")     head(my.mat)          #Ensure that input data is Z-transformed      pca.object <- prcomp(my.mat,center=TRUE,scale.=TRUE)     summary(pca.object)     par(mfrow=c(1,2))     plot(pca.object)     biplot(pca.object)          #The Rotation     pca.object$rotation


title: What is a fuzzy neuron?

text: The only thing I understood about fuzzy neurons is that the neuron's activation function is replaced with some operation used in fuzzy. Other than that, I didn't understand much from science paper I found. Also, it is hard to just Google it because of it's name, which results in listing some fuzzy or ANN or hybrid systems, which is not what I'm looking for. Is there anything more to fuzzy neurons than I have described above? Also, I'm curious, **when and how to use it**? And is there **any other change in neural network structure or calculation that needs implementing** if I'm going to use fuzzy neuron. If someone has a good material on this subject, please link. I really have had a hard time finding anything on this subject. EDIT: Title of the paper I've mentioned: Combining neural networks and fuzzy logic for applications in char.recognition. ByAnne Magaly de Paula Canuto. University of Kent. Can't find link to it anymore, i found it last year. But this is the title. Fuzzy neuron is mentioned in section 1.4.


title: Permutation tests of correlation matrix in R , which correlations are significant?

text: As title, I have a correlation matrix available with cor() function                   corMatrix<-cor(mydata)         corMatrix is a 98 by 98 matrix storing positive and negative correlation values.      I want to assess which correlations values are significant among all the correlations. Permutation test has been suggested to do this. How to perform a permutation test on the correlation matrix ? How to obtain p-values for each correlation value estimated ? Any package suggested ? Or answer with working code will be well appreciated.


title: Assumptions of linear regression

text: One particular assumption is that $e_{i} \sim N(0,\sigma^{2})$ but I was wondering why we use a QQ plot to test this assumption since we have $n$ iid distributions to look at. When we use a QQ plot, aren't we looking for the normality of all the residuals as opposed to one particular $e_{i}$? I am confused about this assumption. Does the QQ plot somehow imply that each $e_{i}$ is normally distributed?


title: Probability distribution satisfying constraints?

text: Reposted from Math.SE: Continuing from this question. Given two random variables $X$ and $Y$ where $X \sim \operatorname{Beta}(a, b)$ and $Y \sim \operatorname{Beta}(c, d)$, I'm looking for a random variable $Z$ with a distribution supported on $[0, 1]$ that satisfies the following constraints:   1. The pdf of $Z$ has reflection symmetry: $f(x;a, b, c, d)=f(1-x;d, c, b, a)$   2. If $\mathrm{E}[X] > \frac{1}{2}$, then $\mathrm{E}[Z] > \mathrm{E}[Y]$   3. If $\mathrm{E}[X] = \frac{1}{2}$, then $\mathrm{E}[Z] = \mathrm{E}[Y]$   4. If $\mathrm{E}[X] < \frac{1}{2}$, then $\mathrm{E}[Z] < \mathrm{E}[Y]$   5. Increasing the expectation of $X$ or $Y$ must never decrease the expectation of $Z$. Is an answer to this question even possible? The answer to my first question does not seem to solve the problem (the expectations are OK, but the support constraint is violated).


title: Which model should I use? Logistic regression?

text: Here is the data I have:   * Response variable : It contains proportions and it takes discrete values 0, 0.2, 0.4, 0.6, 0.8, 1. But there are 109 possible discrete values   * Predictor variable.1: Discrete and ordinal. It contains these values 10, 20, 30, 40.   * Predictor variable.2: Discrete and non-ordinal. It contains these values 'a', 'b', 'c'. Neither normality (checked with Kolmogorov-Smirnov and by looking at a qqplot) nor homoscedasticity (checked with Fligner and by looking at a plot) are respected. Which model should I use in order to infer whether any of the two predictor variable influence my response variable? What about a logistic regression? Would it work?


title: Is it ok to determine early stopping using the validation set in 10-fold cross-validation?

text: I am working on a machine learning experiment comparing the use of multiple different neural network classifiers by applying them on a large number of datasets, using stratified 10-fold cross-validation. I measure the performance as the average of the errors on the validation set (sometimes referred to as test set) of the 10-fold cross-validation procedure. My question is, would it be ok to use this same validation set to do an early stopping of the training procedure? This early stopping would be performed by applying the trained model after each epoch to the validation set and measuring the performance, and if it declines for a number of successive learning epochs, the learning would be halted and we would take the epoch that produced the last good performance. This would be applied to all the different techniques, and across all the different datasets. Is this ok? Or is it statistically inaccurate?


title: Coefficient of variation to compare running time of two algorithms

text: I have developed two new sorting techniques using C language. I need to compare the performance of the both sorting techniques, to see which one is better than another. To do this, I used different input data, from n = 500 to n = 2500. For each data that need to be sorted, I run them 10 times and get the mean, standard deviation as well as coefficient of variation. An example of results are shown below:               n          min          max          mean         SD         CV     500        24.07        24.52        21.28        0.11       0.47     1000       52.41        52.83        52.67        0.13       0.25      Just to ask, can I only run each experiment 10 times or should I run more? If I run them 50 times, the results are very close and I obtain almost the same CV. The above results obtained when I run the first sorting technique. Can I say that 10 times is OK to get the results? I don't know how many samples should I use to test a particular data. Can anyone suggest any good reading books for this particular analysis?


title: Adding weights for highly skewed data sets in logistic regression

text: I am using a standard version of logistic regression to fit my input variables to binary output variables. However in my problem, the negative outputs (0s) far outnumber the positive outputs (1s). The ratio is 20:1. So when I train a classifier, it seems that even features that strongly suggest the possibility of a positive output still have very low (highly negative) values for their corresponding parameters. It seems to me that this happens because there are just too many negative examples pulling the parameters in their direction. So I am wondering if I can add weights (say using 20 instead of 1) for the positive examples. Is this likely to benefit at all? And if so, how should I add the weights (in the equations below). The cost function looks like the following: $$J = (-1 / m) \cdot\sum_{i=1}^{m} y\cdot\log(h(x\cdot\theta)) + (1-y)(1 - \log(h(x\cdot\theta)))$$ The gradient of this cost function (wrt $\theta$) is: $$\mathrm{grad} = ((h(x\cdot\theta) - y)' \cdot X)'$$ Here $m$ = number of test cases, $x$ = feature matrix, $y$ = output vector, $h$=sigmoid function, $\theta$ = parameters we are trying to learn. Finally I run the gradient descent to find the lowest $J$ possible. The implementation seems to run correctly.


title: Item response theory vs. summative Likert scale

text: 1. If the items to be "summed" or combined to create an overall index are collectively the underlying construct (i.e. I am trying to measure compliance to an intervention which has different components), wouldn't combining all the components for IRT violate the assumption of unidimensionality since they represent the intervention itself? If so, do you have a reference for this?   2. Are there any good references that say if Likert scores are summed into a scale and the Cronbach's alpha is relatively high (~0.8), this is appropriate to use as an index?   3. Are there any other tests of validity I should do other than Cronbach's alpha if all I have done is sum the responses? Thanks in advance, any advice is appreciated!


title: Pros and Cons: Methods for Detrending Time Series Data

text: My memory is fuzzy on the advantages and disadvantages of various methods for detrending time-series data. I'm looking for a succinct summary of why and when one should or should not use the following:   * Differenced data   * Log-differenced data   * Error term, after regressing on _only_ a linear or polynomial time series (e.g., 0,1,2,3,...,t)


title: Dealing with few instances in support vector regression

text: What is the best way of going about dealing with few instances in support vector regression, e.g. only approximately 40? Also - is there an optimal way of dealing with outliers in this case of few instances?


title: Does random effect always produce the same result as fixed effect?

text: I have two models - one is including a categorical covariate as a fixed effect, the other includes it as a random effect:               require(nlme)     set.seed(123)     n <- 100     k <- 5     cat <- as.factor(rep(1:k, n))     cat_i <- 1:k # intercept per kategorie     x <- rep(1:n, each = k)     sigma <- 0.2     alpha <- 0.001     y <- cat_i[cat] + alpha * x + rnorm(n*k, 0, sigma)     plot(x, y)          m2 <- lm(y ~ cat + x)     summary(m2)          m3 <- lme(y ~ x, random = ~ 1|cat, na.action = na.omit)     summary(m3)      As you can see, both models `m2` and `m3` produce exactly the same coefficient estimate for x (including SE). Also the residual standard error is the same. The same result is produced when I simulate some missing data:               # simulate missing data     y[c(1:(n/2), (n*k-n/2):(n*k))] <- NA          m2 <- lm(y ~ cat + x)     summary(m2)          m3 <- lme(y ~ x, random = ~ 1|cat, na.action = na.omit)     summary(m3)      So can we say in general that adding effect as random will have the same impact on the other coefficients and the overall inference as adding it as fixed? If not, can you please provide a simple example (or change the provided one) when this fails?


title: Adding random effect influences coefficient estimates

text: I have always been taught that random effects only influence the variance (error), and that fixed effects only influence the mean. But I have found an example where random effects influence also the mean - the coefficient estimate:               require(nlme)     set.seed(128)     n <- 100     k <- 5     cat <- as.factor(rep(1:k, each = n))     cat_i <- 1:k # intercept per kategorie     x <- rep(1:n, k)     sigma <- 0.2     alpha <- 0.001     y <- cat_i[cat] + alpha * x + rnorm(n*k, 0, sigma)     plot(x, y)          # simulate missing data     y[c(1:(n/2), (n*k-n/2):(n*k))] <- NA          m1 <- lm(y ~ x)     summary(m1)          m2 <- lm(y ~ cat + x)     summary(m2)          m3 <- lme(y ~ x, random = ~ 1|cat, na.action = na.omit)     summary(m3)      You can see that the estimated coefficient for `x` from model `m1` is -0.013780, while from model `m3` it is 0.0011713 - both significantly different from zero. Note that when I remove the line simulating missing data, the results are the same (it is full matrix). Why is that? PS: please note I am not a professional statistician, so if you are about to respond with a lot of math then please make also some simple summary for dummies :-)


title: Cumulative predicted probabilities in logistic regression

text: I'm analyzing users' in-game data in order to model whether they're going to be paid user or not. Here's my model:               Logistic Regression Model          lrm(formula = becomePaid ~ x1 + x2 +          x3 + x4 + x5 + x6, data = sn, x = TRUE,          y = TRUE)                                Model Likelihood     Discrimination    Rank Discrim.                                  Ratio Test            Indexes          Indexes            Obs         1e+05    LR chi2    1488.63    R2       0.147    C       0.774          0          99065    d.f.             6    g        1.141    Dxy     0.547          1            935    Pr(> chi2) <0.0001    gr       3.130    gamma   0.586         max |deriv| 8e-09                          gp       0.011    tau-a   0.010                                                    Brier    0.009                                              Coef    S.E.   Wald Z Pr(>|Z|)     Intercept      -6.7910 0.0938 -72.36 <0.0001      x1              0.0756 0.0193   3.92 <0.0001      x2              0.0698 0.0091   7.64 <0.0001      x3              0.0020 0.0002  11.05 <0.0001      x4              0.0172 0.0057   3.03 0.0024       x5              0.0304 0.0045   6.82 <0.0001      x6             -0.0132 0.0042  -3.17 0.0015        And in my model, I created couple of use cases such as:                   test1   test2            x1  8           9     x2  10          10     x3  250        250     x4  6           6     x5  2           2     x6  0           1      Then the probability of user test1 is to turn out to be a paid user is %.07 and % 0.84 for test2. However I want to calculate the cumulative probabilities such as users whose' x1 values are greater than 8, x2 values are between 10 and 20 and so on. Is there any way to calculate this ? Thanks !


title: What do we call multiple testing?

text: I can think of different “types” of multiple testing when using linear models for example:   1. Multiple inferences because we have several dependent variables   2. Multiple inferences because we have several independent variables   3. Looking at the data without making any test. Running a test on only the comparisons that might possibly yield to a significant p.value.   4. Running multiple different tests on the same data. (try a LM, if it is not significant, try a GLM, if it is still not significant try a beta regression, etc.) Wikipedia says: > [...] multiple testing problem occurs when one considers a set of > statistical inferences simultaneously or _infers a subset of parameters > selected based on the observed values_. Does the first part of what wikipedia says encompass my first two points and the part that I put in italics (after the "or") is equivalent to my third point ? Is it correct that my point 4 has nothing to do with what we call multiple testing ? If my question is too blurry, I might rephrase it this way: **When does an issue of multiple testing occur ? How would you categorize (if needed) the possible events of multiple testing ?**


title: Basic proof of mixture models

text: I need to proof that $X$ follows a distribution $F$ with probability $1-p$ and a distribution $G$ with probability $p$ if, and only if, its distribution function is: $(1-p)F + pG$ Can anyone give me some hints about how to proof this? **EDIT2** : this is my _second_ take based on comments: $$F_p(x)=Pr(X\le x)=\int_{-\infty}^x dF_p(x)=(1-p)\int_{-\infty}^xdF(x)+p \int_{-\infty}^xdG(x) = (1-p)F+pG$$ but I am not sure if this is a valid proof


title: What is the meaning of using distributions for statistical inference?

text: In statistics people gather a lot of data (ie heights of people or gene expression levels) to get some insight. Then in order to perform statistical analyses they try to fit their data to a theoretical distribution (ie Normal Distribution) by computing some parameters. How do we know that our data follows such a distribution? If we were able to measure all objects in a population and draw a distribution we could probably get something different, a different shape than the theoretical distribution we thought of. Aren’t our calculations wrong then if we use a theoretical distribution? Please correct me if I am wrong and tell me what you think of this. I think this is a very basic concept in statistics and I have to clarify it.


title: Do I need p-value adjustment for multiple ANOVAs on the same dataset?

text: I have a 5 x 2 design where one of the levels of the first factor is a control, all others are experimental conditions. I'm interested in the interaction between the two factors and especially, if the interaction is present for each of the 4 experimental conditions with the control condition. I conduct 4 separate 2 x 2 ANOVAs where I pair each experimental condition with the control condition in the first factor. This seems to call for p-value adjustment to avoid the multiple-testing problem but what do I have to adjust? Are all p-values adjusted (the two main effects and the interaction)? Or do I only adjust the interaction p-values? Or separately for the main effects? Thanks a lot


title: Analysis of complex sample

text: The company where I work did a survey with a complex sample. Originally, the sample was a 2-step stratified, but due to some problems, we lost one of them (so, let's consider that the sample is stratified only on 1 variable). The biggest part of the questions is like "What you think the company can do in respect of X?" and then there is a list of options (some questions allow only 1 answer, some more than 1, others are Likert-type response options). My question is: Wow can I analyze this data? It would be nice if there were a test like Chi-square for complex sample (probably there exists one, but I don't know) Any software is welcome, but I prefer R and SPSS.


title: Whether to use multivariate F test, df adjusted estimate, or sphericity assumed in repeated measures ANOVA?

text: I have data that involves 2 groups (equal sample size in each) and data for each group over 3 time points (they are actually 3 different monetary reward conditions). I want to investigate within group differences. All time/condition points are important, I don't have a "control" time point. Any opinions on the following:   1. If I am comfortable that sphericity is assumed by Mauchly's Test not being violated, would a multivariate Lambda F test statistic be more appropriate or a sphericity assumed estimate (within-condition estimate from SPSS)?   2. If I choose to go with the understanding that despite a test result saying it is OK, assuming sphericity may still be over-confident, any thoughts on potentially reporting ALL Greenhouse-Geiser estimates regardless of Mauchly's or Lambda? I have read that this may reduce the chance of a Type-1 error without having to assume sphericity or equality of covariance matrices. Perhaps too overcautious? Or, is potentially adjusting df's more invasive than assuming sphericity?


title: What differentiates a population risk model versus Individual risk model?

text: I am trying to determine why people ask if a model of a health care risk (for inpatient stay, or acquiring a disease or any outcome) is a population risk model or an individual risk model. This is in a clinical setting where the outcome is some medical event and covariates are personal health history and demographics. The risk here is then the prob(outcome) as determined by some binary classifier. I understand there may be a different emphasis on the way quality of the model is assessed. For an individual risk model it is perhaps more important to have a small variance of an individual's score and this measure is not captured by AUC of the model or other traditional accuracy measures. Is that the main reason for the difference? If the model is created with observations of people's health history so rows in model are people, then is this necessarily an individual risk model? If not why not? Thanks


title: Some exercises related to sampling and the normal distribution

text: Following cardiac surgery, patients are encouraged to exercise regularly (assume that regular exercise is defined as exercising on 3 or more days per week). A physician suspects that patients exercise regularly immediately following cardiac surgery but tend to reduce, even stop exercising completely, over time. An investigation is planned to estimate the mean number of weeks that patients exercise regularly following cardiac surgery. Assume that the standard deviation (s.d) in the number of weeks cardiac patients exercise regularly following surgery is 6.3 weeks. (a) If a sample of 40 cardiac patients is followed, and the number of weeks in which each patient exercised regularly is recorded, what is the probability that the sample mean will be no more than 1 week higher than the true mean? Under central limit therom,               standard error = s.d./sqrt(n) = 6.3/sqrt(40) = 0.99     X = true mean + 1     Z = (X - true mean) / 0.99 = (true mean + 1 = true mean) / 0.99 = 1/0.99      (b) Find the probability that the sample mean is at least two weeks less than the true mean for a sample of 40 cardiac patients.               X = true mean - 2     Z = (X - true mean) / 0.99 = (true mean -2 - true mean) / 0.99 = -2/0.99      (c) If the sample is increased to 100 cardiac patients, what is the probability that the sample mean will be no more than 1 week higher than the true mean? same as question Can anyone help me check my answers to see whether is correct or not? Please correct and explain if wrong...


title: Use law of total variance to find unconditional variance of overdispersed Poisson?

text: First, I need to prove that the distribution of a RV X, where X|lambda ~ Pois(lambda), and lambda ~ gamma(a, B), is a negative binomial. I know that it _is_ , but why negative binomial instead of another 2 parameter distribution? How do I prove negative binomial is the best alternative in the circumstance of an overdispersed Poisson, basically? Then, I have to calculate V(X). I'm wondering if V(X) will help prove that negative binomial is best to use as an alternative in an overdispersed Poisson?? I looked at this post, and I desperately want to understand @probabilityislogic's math. I am new to notation associated with probability concepts (not _that_ new, though), and I don't know how he/she jumps from line 2 to 3: > $λ_i∼Gamma(α,β)$   >  On doing the integration/mixing over λi, you have:   >  $Y_i(t_i)|αβ∼NegBin(α,p_i)$ where $p_i=\frac{t_i}{t_i+β}$


title: Interpretation of garch results

text: How can I interpret McLeod-Li test results below               > McLeod.Li.test(y = data, gof.lag = 12, plot = FALSE)$p.values      [1] 7.602866e-09 0.000000e+00 0.000000e+00 0.000000e+00 0.000000e+00 0.000000e+00      [7] 0.000000e+00 0.000000e+00 0.000000e+00 0.000000e+00 0.000000e+00 0.000000e+00      And also fitting result .In fitting results, sometimes it says false convergence and sometimes relative function convergence what does it mean?.                ***** ESTIMATION WITH ANALYTICAL GRADIENT *****                     I     INITIAL X(I)        D(I)               1     2.654366e+00     1.000e+00          2     5.000000e-02     1.000e+00          3     5.000000e-02     1.000e+00              IT   NF      F         RELDF    PRELDF    RELDX   STPPAR   D*STEP   NPRELDF          0    1  2.042e+02          1    4  2.041e+02  7.06e-05  2.05e-04  3.0e-03  1.3e+02  1.8e-02  1.30e-02          2    5  2.041e+02  8.08e-05  1.08e-04  2.9e-03  2.0e+00  1.8e-02  6.01e-03          3    6  2.041e+02  1.03e-04  1.03e-04  2.4e-03  2.0e+00  1.8e-02  4.82e-03          4   10  2.037e+02  1.89e-03  1.87e-03  7.1e-02  1.7e-01  3.7e-01  4.86e-03          5   12  2.029e+02  3.77e-03  4.42e-03  1.9e-01  1.9e+00  7.4e-01  5.03e-01          6   15  2.028e+02  6.83e-04  1.23e-03  6.3e-03  4.0e+00  2.6e-02  4.34e-02          7   16  2.028e+02  1.07e-04  1.29e-04  7.5e-03  2.0e+00  2.6e-02  9.17e-03          8   17  2.028e+02  1.47e-04  2.09e-04  1.7e-02  1.9e+00  5.2e-02  5.64e-03          9   18  2.027e+02  1.32e-04  1.68e-04  1.7e-02  1.8e+00  5.2e-02  1.41e-03         10   20  2.027e+02  1.59e-04  2.20e-04  4.3e-02  7.8e-01  1.2e-01  3.48e-04         11   21  2.027e+02  4.71e-05  6.10e-05  2.9e-02  0.0e+00  8.2e-02  6.10e-05         12   22  2.027e+02  1.29e-06  1.37e-06  5.2e-04  0.0e+00  1.8e-03  1.37e-06         13   23  2.027e+02  1.92e-07  5.61e-08  3.6e-04  0.0e+00  9.3e-04  5.61e-08         14   35  2.027e+02 -9.82e-16  1.27e-19  3.3e-15  8.0e+08  9.1e-15  2.68e-10           ***** FALSE CONVERGENCE *****           FUNCTION     2.026871e+02   RELDX        3.290e-15      FUNC. EVALS      35         GRAD. EVALS      14      PRELDF       1.274e-19      NPRELDF      2.677e-10               I      FINAL X(I)        D(I)          G(I)               1    1.240313e+00     1.000e+00     1.250e-03          2    1.234210e-01     1.000e+00     1.351e-05          3    4.425966e-01     1.000e+00     2.553e-03       Another results               ***** ESTIMATION WITH ANALYTICAL GRADIENT *****                     I     INITIAL X(I)        D(I)               1     3.736575e+00     1.000e+00          2     5.000000e-02     1.000e+00              IT   NF      F         RELDF    PRELDF    RELDX   STPPAR   D*STEP   NPRELDF          0    1  2.061e+02          1    2  1.764e+02  1.44e-01  1.61e+00  1.3e-01  3.3e+02  1.0e+00  2.69e+02          2    5  1.229e+02  3.03e-01  1.99e-01  3.7e-01  7.8e-01  2.0e+00  1.12e+00          3    7  1.070e+02  1.30e-01  9.76e-02  1.3e-01  3.7e+00  4.0e-01  2.19e+02          4    9  6.528e+01  3.90e-01  3.24e-01  4.3e-01  8.1e+03  8.0e-01  5.81e+03          5   11  5.613e+01  1.40e-01  1.37e-01  1.1e-01  8.8e+00  1.6e-01  3.40e+05          6   13  5.456e+01  2.79e-02  2.88e-02  2.1e-02  3.2e+01  3.2e-02  4.02e+02          7   15  5.217e+01  4.39e-02  4.78e-02  4.3e-02  2.6e+00  6.4e-02  1.08e+02          8   17  5.161e+01  1.07e-02  1.94e-02  3.7e-02  2.6e+00  5.7e-02  6.86e+01          9   18  5.132e+01  5.52e-03  8.37e-03  3.4e-02  2.0e+00  5.7e-02  4.83e+01         10   19  5.129e+01  7.58e-04  4.77e-03  3.0e-02  2.0e+00  5.7e-02  1.14e+01         11   21  5.123e+01  1.02e-03  2.86e-03  1.1e-02  2.0e+00  2.7e-02  5.80e-01         12   22  5.121e+01  3.58e-04  5.08e-04  1.5e-02  2.0e+00  2.7e-02  2.63e-03         13   25  5.121e+01  3.94e-05  6.44e-05  1.3e-03  3.8e+00  2.4e-03  4.69e-02         14   26  5.121e+01  2.55e-05  2.88e-05  1.3e-03  3.1e+00  2.4e-03  1.50e-02         15   27  5.121e+01  2.39e-05  5.20e-05  2.3e-03  2.0e+00  4.8e-03  8.22e-03         16   28  5.121e+01  5.18e-05  1.27e-04  5.4e-03  1.9e+00  9.6e-03  2.78e-03         17   30  5.121e+01  6.65e-06  2.50e-05  1.2e-03  1.8e+00  3.0e-03  8.02e-05         18   31  5.121e+01  4.39e-06  6.27e-06  1.7e-03  1.3e+00  3.0e-03  6.91e-06         19   33  5.121e+01  5.81e-07  2.26e-06  5.6e-04  1.5e+00  1.1e-03  4.17e-06         20   34  5.121e+01  5.53e-07  5.82e-07  2.0e-04  0.0e+00  4.6e-04  5.82e-07         21   35  5.121e+01  3.22e-09  3.35e-09  4.8e-05  0.0e+00  8.6e-05  3.35e-09         22   36  5.121e+01  8.01e-12  8.10e-12  2.4e-06  0.0e+00  4.2e-06  8.10e-12           ***** RELATIVE FUNCTION CONVERGENCE *****           FUNCTION     5.120700e+01   RELDX        2.408e-06      FUNC. EVALS      36         GRAD. EVALS      23      PRELDF       8.103e-12      NPRELDF      8.103e-12               I      FINAL X(I)        D(I)          G(I)               1    2.113519e-01     1.000e+00     1.349e-05          2    8.757729e-01     1.000e+00     2.232e-06


title: Fitting a distribution for frequency

text: I have a list of objects and their frequencies of occurrence in R                   Names <- c("a","b","c","d","e")         Freq <- c(12,45,67,100,1 ...)      I want to rank the prominence of these objects based on the frequency. Is there a way to fit a distribution on the frequencies and weight them to get ranks instead of simply sorting them in descending order?


title: Using MLE in R stats4

text: I have been trying to estimate the MLE for my joint posterior. I'm using R and the package stats4. I have 14 parameters and two of them are $\geq 0$, which I did not know how to implement (and I was creating NaN due to the minus log posterior required in for the mle function) and I just made it return very high value (1000) if either of the parameters were negative. Is this the right way to solve this problem? As I was forced to change my prior each time (because MLE told me that my prior estimates were way to high) and I find these nonegative parameters going down to were low numbers (0.001 and 0.01) which did not seem right and at each iteration way below my repeatly suggested prior. Also, since I didn't have the exact posterior due to the structure of the model and I tried to scale it such that the point estimate from the mle function plugged in the log joint posterior had the value 0. Is this approximation okay for this function? Regards, Raxel.


Let $X_1, X_2, ..., X_n$ be a random iid sample from a population with mean $\theta$. Now I am wondering about the intuition behind $E(X_1| \overline X ) = \overline X$, the sample mean. If we just consider $X_1$ (or any $X_i$ for that matter) we have that $E(X_1) = \theta$ as the expected value of any random observation from a population will be the population mean. Now given that we know $\overline X$ how come that changes what we expect to get for $X_1$? $X_1$ is still the same random observation from the population as before...but it seems that knowing the sample mean 'overrides' what we expect to get for an observation...is this correct?


title: Parameter estimation based on MLE estimate of another parameter

text: There is a situation explained below where I intend to apply MLE. The problem statement is that I am estimating a measure $X$. This measure is obtained my Maximum Likelihood estimation technique. I will apply this measure to an AR(2) model to find the parameters of AR model. The way the measure is applied and calculated is explained briefly. Considering a linear auto regressive (AR) model whose parameters are needed to be found that is corrupted by $\eta(t)$ white noise. Question : What is the effect of noise on MLE?


title: Loss function for linear regression with calculus of variations

text: I'm struggling with mathematics behind linear regression. In the following lines I pasted the text from the book Pattern Recognition and Machine Learning (p. 46) where author derives the regression function $\mathbb{E}_{t} [t | \mathbf{x}]$. I want to understand the procedure from the equation (2) to the final result. Could somebody please provide me some useful pointers (and/or links) which concept from the calculus of variations should I study. The average, expected, loss is given by $$ \mathbb{E}[L] = \int \int L(t, x (\mathbf{x})) p (\mathbf{x}, t) \, d\mathbf{x} \, dt. \tag{1} $$ A common choice of loss function in linear regression is the squared loos given by $L (t, y(\mathbf{x})) = \{ y (\mathbf{x}) - t \}^{2}$. In this case, the expected loss can be written as $$ \mathbb{E}[L] = \int \int \{ y (\mathbf{x}) - t \}^{2} p (\mathbf{x}, t) \, d\mathbf{x} \, dt. \tag{2} $$ Our goal is to choose $y (\mathbf{x})$ so as to minimize $\mathbb{E} [L]$. We can do this using the calculus of variations to give $$ \dfrac{\delta \mathbb{E} [L]}{\delta y (\mathbf{x})} = 2 \int \{ y (\mathbf{x}) - t \} p (\mathbf{x}, t) \, dt = 0. \tag{3} $$ Solving for $y (\mathbf{x})$, and using the sum and product rules of probability, we obtain $$ y (\mathbf{x}) = \dfrac{\int tp (\mathbf{x}, t) \, dt}{p (\mathbf{x})} = \int t p (t | \mathbf{x}) \, dt = \mathbb{E}_{t} [t | \mathbf{x}] \tag{4} $$


title: Normalization against Covariates

text: I have a list of parameters which correlate with 1-2 covariates that I want to control for. Following normalization, I wanted to do comparisons between groups, correlation analysis and probably use some of them as features in a prediction model. The problem is that I don't think I can use the normalized values following ANOVA for regression analysis (and so on), since they will be by default correlated as a result of the ANCOVA. Do you agree (or not)? What else would you suggest?


title: Comparing densities of a feature for different classes when the feature is irrelevant to one class

text: Let us suppose that I have a number of features. I design pdfs for every feature and every class, some of them by smoothing some histogram of training samples, others just by introducing the prior knowledge on how the feature should look like. Now I want to know which is the most likely class for a given observation, based on the pdfs. I can state the following: $$P(p_i\in C_i | F_i) = \displaystyle\frac{P(F_i | p_i\in C_i ) P(p_i\in C_i)}{\displaystyle \sum_{C_j} P(F_i | p_i\in C_j ) P(p_\in C_j)}$$ being $p_i$ a particular observation with feature $F_i$, $C_i$ the class I'm assigning to it. What I've doing so far is just replacing $P(\cdot|\cdot)$ by its density function $\delta$ given that $P(\cdot|\cdot)$ is about $2\epsilon \delta(\cdot|\cdot) $ for a neighborhood $[F_i-\epsilon,F_i+\epsilon]$ and I can cancel out the $2\epsilon$, so I suppose it's just legitimate to directly take the value at the density function. Now, there's one feature in particular in which I want to convey that I just don't care about its value in a very wide range, for one of the classes. It's an area, so say the area can be between zero and a large value. In other words, in this class we can find any area. In another class, there is a shorter dynamic range of the values the feature should take. The problem that I find is that the result of the formula before is always lower for the 'don't care' class because the uniform distribution extends far away and hence the value of the pdf is low at every point. My desired behavior would be that since I don't care much about the area, there should be a more fair competition among the two classes, and probably some other features in the classification process will make the final decision, but this way I just say that the problematic class is very unlikely for any input of this feature. I could tune $P(p\in C_i)$ to be higher for the problematic class. However, in other features this is not what I expect. Any clues on how to handle this or pose this problem?