title: How to quantify intraindividual variability between paired measurements?

text: I've collected data measuring several different parameters in several subjects at two time-points (t1 and t2). X and Y are meausurements of two anatomical structures; X is an area measured in mm2, Y is a volume measured in mm3. For each one I can describe at least 4 variables: measuments at t1 and t2 (eg. X1, X2), absolute difference between the two measurements (X2-X1), percentual difference between the two. I want to compare the "reliability" of X and Y in term of intraindividual variability. My intuitive choice would be to compare standard deviation or coefficient of variation between the two (or their variation). But:   * In the single time point measure (say, X1 vs Y1), a low SD/CV tells me only that different individuals tend to have roughly the same values, while high SD/CV tells me that different individuals tend to have very different values.   * Considering absolute difference between the two measurements does not seems to have much sense.    * Considering percentual difference between the two, low SD/CV in e.g. X tells me only that, in a population, X tend to change of roughly the same % in all the individuals, while an high SD/CV tells me that there is more variability in the entity of the change between individuals. So, I'm not sure which statistic will tell me what I want to know. I was thinking of Pearson's r (or Spearman's rho) to quantify correlation between X1-X2 vs correlation between Y1-Y2 (if there is a stronger correlation between X1 and X2 than with Y1-Y2, it means X1 and X2 tend to increase/decrease in the same direction more thant Y1 and Y2 do, so they is less intraindividual variability).   * Am I correct in this assumption?   * If I were to calculate this in R, is it just a matter of cor.test(X1,X2) cor.test(Y1,Y2) and then compare results?


What is the difference between finite and infinite variance ? My stats knowledge is rather basic; Wikipedia / Google wasn't much help here.


title: MATLAB interperetation of Neural networks

text: I am new to Neural networks and I am trying to build a custom neural network using the NN toolbox in MATLAB.I am using the "create custom neural network function". Now, I find the neural network interpretation of MATLAB very confusing. I am looking to create a network which looks like this - ![NN](http://i.stack.imgur.com/t2wmJ.png) How do I do this?when I try to make Neural networks in MATLAB, the result is somewhat like this - ![NN2](http://i.stack.imgur.com/yEpNo.png) I am having problems interpreting this in terms in neurons in a neural network. I hope the problem is understandable. Thanks in advance.


title: Information Criteria: penalization for number of parameters or for within sample?

text: After I read chapter 7 of the new edition of Bayesian Data Analysis1, I have come to understand that while Information Criteria like DIC and WAIC are:   * A way measuring the adjust of the model and penalize by the effective number of parameters.    They are also:   * A way to get the expected out-of-sample (pointwise) log predictive density by starting with the within-sample (pointwise) log predictive density and correct its bias (with the measure of effective number of parameters).    Although I does seem intuitive that the more parameters we have in the model, the more we tend to overestimate how well the model explains with-sample predictive accuracy, I still can't see why these two things are exactly the same thing, specially how the number of parameters gets us from within-sample to expected out-of-sample. (Or I am wrong, and if so, where/why) [1] Gelman, A.; Carlin, J. B.; Stern, H. S.; Dunson, D. B.; Vehtari, A. & Rubin, D. B. (2014). Bayesian data analysis. CRC Press, third edition.


title: Statistics that shows how two ranges of numbers overlap?

text: I have collected ratings on a scale from 0 to 100 for two different conditions: A and B. I want to show now, that the ratings collected under both conditions are very similar and that their range overlap. How can I do that best? I thought about drawing a box plot for both conditions, but I would prefer a number / numbers to show.


title: Sampling to maximize model accuracy

text: Suppose you have a relatively small random sample and have a corresponding model $\ Y$ ~ $\operatorname{Bernoulli}(p_i) $ $\ \operatorname{logit}( \hat{p_i} )=\hat{\beta}*X$ and now want to draw a sample of a fixed size in order to maximize the accuracy over the entire universe of the $\hat{p_i}$ generated after I refit the model with the new data. The idea being that there might be subgroups that are a large percent of the population but have extremely predictable opinion, meaning that you don't need to sample very many of them to get an accurate model. This, or some variant, has to be a pretty well known problem. Is there any literature I can read about it?


title: JAGS prior-based prediction

text: I'm trying to use rjags to predict what the curve-maximums will be for different groups of incomplete data (one metric x and one metric y). Here is an example of the indata: ![enter image description here](http://i.stack.imgur.com/315aR.png)               in.data <- structure(list(x = c(10.4055555555556, 16.1694444444444, 25.25,                                      33.3305555555556, 37.4736111111111, 31.2805555555556, 24.3347222222222,                                      24.1638888888889, 27.8222222222222, 34.0055555555556, 31.6611111111111,                                      33.8666666666667, 26.6486111111111, 25.3986111111111, 24.8611111111111,                                      15.2125, 6.33194444444444, 11.7541666666667, 15.6694444444444,                                      27.9430555555556, 34.1180555555556, 36.5092592592593, 34.4763888888889,                                      32.8055555555556, 26.9041666666667, 27.9319444444444, 28.7319444444444,                                      30.3638888888889, 31.3861111111111, 26.4722222222222, 27.0152777777778,                                      26.6333333333333, 17.3902777777778, 8.09861111111111, 12.9305555555556,                                      21.8277777777778, 30.8652777777778, 36.9486111111111, 33.5722222222222,                                      29.1652777777778, 25.6444444444444, 28.7222222222222, 32.9319444444444,                                      32.8375, 31.9402777777778, 29.4569444444444, 26.9013888888889,                                      25.7152777777778, 24.9638888888889, 14.3722222222222, 7.36527777777778,                                      13.4722222222222, 20.7263888888889, 31.6388888888889, 41.6555555555556,                                      49.0513888888889, 39.1333333333333, 33.5847222222222, 28.2722222222222,                                      28.2583333333333, 31.2361111111111, 33.1569444444444, 29.9166666666667,                                      27.925, 26.9680555555556, 24.7611111111111, 18.8319444444444,                                      6.33055555555556, 9.87083333333333, 14.5541666666667, 18.5041666666667,                                      23.0652777777778, 21.8875, 20.7291666666667, 14.7375, 11.4083333333333,                                      12.7930555555556, 12.8013888888889, 16.8277777777778, 14.8263888888889,                                      13.3222222222222, 9.74305555555556, 6.42638888888889, 6.275,                                      9.76666666666667, 9.98472222222222, 8.5, 11.4958333333333, 10.9777777777778,                                      13.3083333333333, 8.69722222222222, 6.24305555555556, 5.11111111111111,                                      6.41527777777778, 8.04583333333333, 9.46111111111111, 11.0333333333333,                                      7.3, 9.42916666666667, 8.86805555555556, 9.84305555555556, 10.8555555555556,                                      7.32916666666667, 5.77361111111111, 7.52083333333333, 9.99305555555556,                                      10.9875, 10.0888888888889, 11.1, 13.8277777777778, 15.3458333333333,                                      9.74166666666667, 6.59305555555556, 7.28611111111111, 5.95833333333333,                                      9.20972222222222, 8.85138888888889, 11.6277777777778, 12.9680555555556,                                      12.8541666666667, 14.1041666666667, 11.1236111111111, 6.77777777777778,                                      6.48055555555556, 9.77361111111111, 8.14444444444444, 7.92222222222222,                                      7.85555555555556, 6.63611111111111, 5.40416666666667, 0, 1, 2,                                      3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19,                                      20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35,                                      36, 37, 38, 39, 40, 41, 42, 43, 44, 45, 46, 47, 48, 49, 50, 51,                                      52, 53, 54, 55, 56, 57, 58, 59, 60, 61, 62, 63, 64, 65, 66, 67,                                      68, 69, 70, 71, 72, 73, 74, 75, 76, 77, 78, 79, 80, 81, 82, 83,                                      84, 85, 86, 87, 88, 89, 90, 91, 92, 93, 94, 95, 96, 97, 98, 99,                                      100), y = c(1.36111111111111, 2.08805555555556, 2.83666666666667,                                                  3.17583333333333, 3.53194444444444, 1.4125, 2.44166666666667,                                                  2.45777777777778, 2.78527777777778, 3.17611111111111, 3.02138888888889,                                                  2.60861111111111, 2.88805555555556, 2.79527777777778, 2.575,                                                  1.77333333333333, 0.937777777777778, 1.57027777777778, 1.69888888888889,                                                  2.74, 3.32888888888889, 2.6975, 3.42444444444444, 3.40555555555556,                                                  2.77277777777778, 2.6425, 2.96166666666667, 3.25472222222222,                                                  3.23777777777778, 2.70472222222222, 2.76805555555556, 2.39694444444444,                                                  1.97083333333333, 1.00722222222222, 1.58055555555556, 2.52138888888889,                                                  3.16083333333333, 3.30555555555556, 3.61805555555556, 3.11472222222222,                                                  2.7825, 3.11416666666667, 3.25944444444444, 3.12666666666667,                                                  3.50555555555556, 3.07972222222222, 2.92666666666667, 2.45361111111111,                                                  2.33055555555556, 1.61555555555556, 0.814444444444444, 1.41722222222222,                                                  2.30222222222222, 3.27166666666667, 3.70638888888889, 5.32527777777778,                                                  3.8375, 3.42, 2.985, 2.99194444444444, 3.355, 3.36472222222222,                                                  3.07388888888889, 3.095, 3.00222222222222, 2.61583333333333,                                                  2.38277777777778, 0.776111111111111, 1.15611111111111, 1.53361111111111,                                                  2.16611111111111, 2.62861111111111, 2.60777777777778, 2.46111111111111,                                                  1.91111111111111, 1.35083333333333, 1.50833333333333, 1.50722222222222,                                                  1.63555555555556, 1.61527777777778, 1.53333333333333, 1.09055555555556,                                                  0.396666666666667, 0.401666666666667, 0.547777777777778, 0.594722222222222,                                                  0.356388888888889, 0.790277777777778, 0.510555555555556, 0.823611111111111,                                                  0.563611111111111, 0.260277777777778, 0.146944444444444, 0.363333333333333,                                                  0.369444444444444, 0.585, 0.740833333333333, 0.320833333333333,                                                  0.765277777777778, 0.698611111111111, 0.6625, 0.764722222222222,                                                  0.611388888888889, 0.276388888888889, 0.518611111111111, 0.735,                                                  0.779722222222222, 0.732222222222222, 0.882222222222222, 1.02444444444444,                                                  1.06638888888889, 0.730555555555556, 0.469444444444444, 0.565277777777778,                                                  0.573333333333333, 0.720277777777778, 0.780277777777778, 0.898888888888889,                                                  1.02305555555556, 0.910277777777778, 1.01194444444444, 0.884722222222222,                                                  0.541388888888889, 0.469166666666667, 0.709166666666667, 0.599166666666667,                                                  0.502777777777778, 0.523333333333333, 0.547222222222222, 0.243333333333333,                                                  0.0156494522691706, 0.264475743348983, 0.510172143974961, 0.752738654147105,                                                  0.992175273865415, 1.22848200312989, 1.46165884194053, 1.69170579029734,                                                  1.91862284820031, 2.14241001564945, 2.36306729264476, 2.58059467918623,                                                  2.79499217527387, 3.00625978090767, 3.21439749608764, 3.41940532081377,                                                  3.62128325508607, 3.82003129890454, 4.01564945226917, 4.20813771517997,                                                  4.39749608763693, 4.58372456964006, 4.76682316118936, 4.94679186228482,                                                  5.12363067292645, 5.29733959311424, 5.4679186228482, 5.63536776212833,                                                  5.79968701095462, 5.96087636932707, 6.1189358372457, 6.27386541471049,                                                  6.42566510172144, 6.57433489827856, 6.71987480438185, 6.8622848200313,                                                  7.00156494522692, 7.1377151799687, 7.27073552425665, 7.40062597809077,                                                  7.52738654147105, 7.6510172143975, 7.77151799687011, 7.88888888888889,                                                  8.00312989045383, 8.11424100156495, 8.22222222222222, 8.32707355242566,                                                  8.42879499217527, 8.52738654147105, 8.62284820031299, 8.7151799687011,                                                  8.80438184663537, 8.89045383411581, 8.97339593114241, 9.05320813771518,                                                  9.12989045383412, 9.20344287949922, 9.27386541471049, 9.34115805946792,                                                  9.40532081377152, 9.46635367762128, 9.52425665101721, 9.57902973395931,                                                  9.63067292644758, 9.679186228482, 9.7245696400626, 9.76682316118936,                                                  9.80594679186228, 9.84194053208138, 9.87480438184664, 9.90453834115806,                                                  9.93114241001565, 9.95461658841941, 9.97496087636933, 9.99217527386542,                                                  10.0062597809077, 10.0172143974961, 10.0250391236307, 10.0297339593114,                                                  10.0312989045383, 10.0297339593114, 10.0250391236307, 10.0172143974961,                                                  10.0062597809077, 9.99217527386542, 9.97496087636933, 9.95461658841941,                                                  9.93114241001565, 9.90453834115806, 9.87480438184664, 9.84194053208138,                                                  9.80594679186228, 9.76682316118936, 9.7245696400626, 9.679186228482,                                                  9.63067292644758, 9.57902973395931, 9.52425665101721, 9.46635367762128,                                                  9.40532081377152), Group_ID = c("Group_1", "Group_1", "Group_1",                                                                                "Group_1", "Group_1", "Group_1", "Group_1", "Group_1", "Group_1",                                                                                "Group_1", "Group_1", "Group_1", "Group_1", "Group_1", "Group_1",                                                                                "Group_1", "Group_1", "Group_1", "Group_1", "Group_1", "Group_1",                                                                                "Group_1", "Group_1", "Group_1", "Group_1", "Group_1", "Group_1",                                                                                "Group_1", "Group_1", "Group_1", "Group_1", "Group_1", "Group_1",                                                                                "Group_1", "Group_1", "Group_1", "Group_1", "Group_1", "Group_1",                                                                                "Group_1", "Group_1", "Group_1", "Group_1", "Group_1", "Group_1",                                                                                "Group_1", "Group_1", "Group_1", "Group_1", "Group_1", "Group_1",                                                                                "Group_1", "Group_1", "Group_1", "Group_1", "Group_1", "Group_1",                                                                                "Group_1", "Group_1", "Group_1", "Group_1", "Group_1", "Group_1",                                                                                "Group_1", "Group_1", "Group_1", "Group_1", "Group_1", "Group_1",                                                                                "Group_1", "Group_1", "Group_1", "Group_1", "Group_1", "Group_1",                                                                                "Group_1", "Group_1", "Group_1", "Group_1", "Group_1", "Group_1",                                                                                "Group_1", "Group_2", "Group_2", "Group_2", "Group_2", "Group_2",                                                                                "Group_2", "Group_2", "Group_2", "Group_2", "Group_2", "Group_2",                                                                                "Group_2", "Group_2", "Group_2", "Group_2", "Group_2", "Group_2",                                                                                "Group_2", "Group_2", "Group_2", "Group_2", "Group_2", "Group_2",                                                                                "Group_2", "Group_2", "Group_2", "Group_2", "Group_2", "Group_2",                                                                                "Group_2", "Group_2", "Group_2", "Group_2", "Group_2", "Group_2",                                                                                "Group_2", "Group_2", "Group_2", "Group_2", "Group_2", "Group_2",                                                                                "Group_2", "Group_2", "Group_2", "Group_2", "Group_2", "Group_2",                                                                                "Group_2", "Prior", "Prior", "Prior", "Prior", "Prior", "Prior",                                                                                "Prior", "Prior", "Prior", "Prior", "Prior", "Prior", "Prior",                                                                                "Prior", "Prior", "Prior", "Prior", "Prior", "Prior", "Prior",                                                                                "Prior", "Prior", "Prior", "Prior", "Prior", "Prior", "Prior",                                                                                "Prior", "Prior", "Prior", "Prior", "Prior", "Prior", "Prior",                                                                                "Prior", "Prior", "Prior", "Prior", "Prior", "Prior", "Prior",                                                                                "Prior", "Prior", "Prior", "Prior", "Prior", "Prior", "Prior",                                                                                "Prior", "Prior", "Prior", "Prior", "Prior", "Prior", "Prior",                                                                                "Prior", "Prior", "Prior", "Prior", "Prior", "Prior", "Prior",                                                                                "Prior", "Prior", "Prior", "Prior", "Prior", "Prior", "Prior",                                                                                "Prior", "Prior", "Prior", "Prior", "Prior", "Prior", "Prior",                                                                                "Prior", "Prior", "Prior", "Prior", "Prior", "Prior", "Prior",                                                                                "Prior", "Prior", "Prior", "Prior", "Prior", "Prior", "Prior",                                                                                "Prior", "Prior", "Prior", "Prior", "Prior", "Prior", "Prior",                                                                                "Prior", "Prior", "Prior", "Prior")), .Names = c("x", "y", "Group_ID"                                                                                ), row.names = c(NA, -231L), class = "data.frame")          library(ggplot2)     p <- ggplot(in.data, aes(x=x, y=y, group=Group_ID, color=Group_ID))     p <- p + geom_point()     p      Group_1 and Group_2 will follow the same parabolic shape as the prior but, as can be seen in the graph, will have lower maximums when x increase. Can Bayesian Regression be used for this? If so, how could the model for rjags look? The closest I think I've gotten is this, but can't get it to work:               N[i] ~ dpois(z[i])     log(z[i]) < beta0 + beta1 * x[i] + beta2 * x[i]^2      I've studied the truly great book by Kruschke (Doing Bayesian Data Analysis) where chapter 16 explains many steps, but when I want to expand from linear y = b1 + b2 * x to a more advanced curve-shape I get stuck.


title: Comparison of multiple diagnostic tests with paired data

text: Good Day All, I am attempting to compare three diagnostic tests against a gold-standard and for simplicity sake, let us have all the outcomes be binary (disease/no disease). In this study study, all participants received the four tests (e.g. everyone got the gold-standard and three competing tests). I see that McNemar's Test is a great way for comparing two tests (A & B) by analyzing the homogeneity of how often the agree/disagree with the gold- standard. What about three tests against the gold standard ? It is my understanding that Bhapkar's Test and/or Maxwell-Stuart Tests are generalizations of McNemar's Test but for more than two outcomes, not more than two tests. Are there other generalizations of other tests that I can use ? Thanks in advance !


title: Assessing spurious effect of a third variable on the relationship between a response and a predictor

text: Say, I have a biomarker that is strongly associated to a gene. This biomarker is also strongly associated to another trait, like glucose, but the gene is not.   * If I perform a regression between the glucose and biomarker + gene I get the biomaker and gene both significant: Is this a spurious effect?    * And what if I add an interaction term between the biomarker and gene on the glucose (biomarker*Gene+Gene+biomarker) and all the terms are significant?  What does it mean if when I add a third variable (biomarker) in the regression the second (gene) becomes significant all of the sudden? Does it mean that the second variable is then significantly associated with the dependent?


title: Need help in automatically computing optimal hyper-parameters for a clustering algorithm

text: I have an algorithm that take as input some data (that are continuously arriving) and 3 or 4 parameter values that should be specified by the user. At the and of execution (or periodically during execution) we can evaluate the performance to see if the chosen parameter values were good or not, but the evaluation function is expensive (time and space consuming). There is no way to know _a priori_ what are the best matching parameter values, unless we perform many execution and evaluations until we find the parameter values that gives a good results. Note that even if we try to variate the parameter values to find the best one by executing the algorithm many times on the same data, these values might be not good for another set of data (this is one reason of why I need to adapt them automatically during execution). EDIT1: Some additional information: the algorithm is supposed to do **online clustering on a continuously arriving data stream.** The parameters are mainly some threshold values or used to compute some adaptive threshold values needed by the algorithm. EDIT2: Some additional information regarding the desired algorithm: Each time a new data x is given to the algorithm, it will see if x should be assigned to an existing cluster (i.e. adapting an existing cluster representative n), or if a new representative based on x should be created. To do that, each representative has a threshold value, if distance between x and nearest representative is less than the threshold of this representative, then x is assigned to this representative, else a new representative will be created based on x. Now some parameters are involved at this level: I assign a default initial threshold value to each new created representative (this is one parameter); when a new data is assigned to an existing representative (cluster) the threshold of this representative is updated, this is what I actually call "adaptive threshold", its defined by "the mean of distances between the representative and its currently associated data" + param * "the standard deviation of this distances" (param here is another parameter). Then I have some other parameters that influence the final number of representatives, like for example the age beyond which an old representative should be deleted ... So here I gave 3 examples of the parameters used by the algorithm. **Questions:**   1. Is there any solution to automatically adapt the parameter values during the execution ?    2. Is it better to try to define a function that adapts these parameters and this function depends on a less sensitive parameter(s) ? Any idea is welcome. Note that my problem is not the clustering algorithm itself; it gives "good results" iif "the good parameter values" are chosen, thus I wanted to know if there anyway to automatically adapt them.


title: Weighted kernel density plot in R

text: I've subsetted and plotted my unweighted data, but I don't see how to make use of my weighting variable "FactBx$expwgt"               #####################################################     plot( density(FactB1$BAV_DIST),       #### xlim=range( c(FactB1$BAV_DIST, FactB2$BAV_DIST, FactB3$BAV_DIST) ),           xlim = c(0, 10),          main = "Density Plots",          xlab = "miles",          col = 2  )          lines(density(FactB2$BAV_DIST), col=3)     lines(density(FactB3$BAV_DIST), col=4)           lines(density(FactB6$BAV_DIST), col=5)            TextVect <- c("Walk","Bicycle","Auto{Dropped}","Auto{Parked}")     ColorVect <- c(2,3,4,5)                  LineVect <- c(0,1,2,3)                    mtext(TextVect, side = 1, line = LineVect, col = ColorVect)     ###########################################################      I've been using the survey package to do weighted summaries, but I don't see how to do weighted density plots.


title: Question on term used for this kind of questionnaire item and type of analysis followed

text: I have a question on the term used for the type of the item that I would like to use and how it is "usually" analysed. Let's assume that we have two items and the one is linked to the other, namely, the one is the answer in a multiple choice question and the second on is quantifies an attribute, for example, "did you find the lesson interesting" then a LIKERT scale is provided and then another item that asks whether finding the lesson interesting was important for the respondent (something like classification/priotirisation). So an item attaches somehow an additional quantity to the first one. These kinds of "pair of items" are used in customer satisfaction research usually... Do you know how they are called? How are they analysed? Would Kullback-Leibler divergence be needed?


title: Is a log transformation a valid technique for t-testing non-normal data?

text: In reviewing a paper, the authors state, "Continuous outcome variables exhibiting a skewed distribution were transformed, using the natural logarithms, before t tests were conducted to satisfy the prerequisite assumptions of normality." Is this an acceptable way to analyze non-normal data, particularly if the underlying distribution is not necessarily lognormal? This may be a very stupid question, but I have not seen this done before....


title: Definition of "Natural Images" in the context of machine learning

text: Whilst reading up on the Deep learning literature, I noticed that a few variations on the standard network structure that were created specifically to better model "Natural/Real Images". For example, this paper says that: > Deep belief nets have been successful in modeling handwritten characters, > but it has proved more difﬁcult to apply them to real images. ... [Our > network] can better model the covariance structure of natural images. The paper seems to imply that only **real or natural** images have a rich local covariance structure. If that's the case, then would a screenshot of a videogame count as a natural image? A digital painting? and in any case in there an algorithmic way to test for this "Naturalness"? This may sound open-ended, so let eplicitly operationalize it: **when given a training set of 2D pixel arrays, how to do determine if you should use a standard network structure or a structure intended for natural/real images?**


title: What are the primary differences between Taxometric analyses (e.g., MAXCOV, MAXEIG) and Latent Class analyses?

text: Recent research has attempted to determine if certain psychological constructs are latently dimensional or taxonic (i.e., including taxons or classes). For example, researchers may be interested in finding out if there is a certain "class" of people who are more likely to develop chronic pain after an injury, or if the risk of developing chronic pain is better conceptualized as dimensional ranging from limited risk to extremely high risk. I've noticed that researchers attempt to answer these types of questions using two types of analyses: Taxometric analyses (MAMBAC, MAXEIG, MAXCOV) typically conducted in R, and Latent Class analyses. **Here are some examples of Taxometric studies:**   * Taxonicity of anxiety sensitivity: a multi-national analysis   * Assessing the latent structure of the intolerance of uncertainty construct: An initial taxometric analysis **Here are some examples using Latent Class analyses:**   * Empirically derived latent classes of tobacco dependence syndromes observed in recent-onset tobacco smokers: epidemiological evidence from a national probability sample survey   * The structure of posttraumatic stress disorder: latent class analysis in 2 community samples **Here are my questions:**   1. In English, what are the primary differences between these two types of analyses? If possible, elaborate if they answer different questions and how they are analytically (mathematically) different.   2. Which one is better for answering the type of question I highlighted in my "introduction", and why? Perhaps this is really unanswerable at this point.  Also please share any information you feel may be relevant to this topic. I have a feeling I will have follow-up questions!


title: What to do when a NN learns a static representation for a large chunk of your dataset?

text: So I'm currently training a stack of RBMs (Restricted Boltzmann Machines), eventually to be build into a DBN (Deep Belief Net), on a set of gray-scale images featuring objects placed in different locations on a complex background. After training for a good long while, I tried generating from the RBMs and found that for a large percentage of the generated images were identical; they all looked like the background. This makes sense: the background is static across the training set, so learning a high fidelity representation of that is a better idea (from the network's perspective) than learning a low fidelity representation of the objects, since they appear in different positions. My question is this: **I want to learn about the whole image, not just the background; what should I change in order get this to happen?** I could always lower the learning rate/momentum and keep training, but I've love some opinions as to whether or not that is the best solution. **Is there anyway to utilize this learned representation of the image background to learn a better representation of the foreground?** I've thought of just subtracting it out, but that seems rather ad hoc. P.S. I've had this problem in the past with traditional autoencoders (i.e. three-layer backpropagation net predicting the input pattern), so I know it isn't unique to deep learning approaches or the particular dataset, just any dataset that contains a static signal amongst more varied information.


title: Is there a preferred panel chart layout?

text: I have a panel chart composed of 4 panels, each with 12 bars that essentially show data from a pivot table. The default layout has the first 11 rows show individual points, with the 12th row showing a total of the entire column. Similarly, the final column shows a total of the 3 preceding columns. The result is the total for all values is bottom right, like in pivot table. My question is whether there's a benefit or preferred method to rearranging the panels and data points so that the aggregate data is top left in the chart, making it the first value seen (in traditional left-right, top-bottom reading), or keeping the more traditional mathematic approach with totals at the bottom and to the right? Thanks for your help.


title: Can we make the Irwin-Hall distribution more general?

text: I need to find a symmetric low-kurtosis distribution class, which includes the uniform, the triangular and the normal Gaussian distribution. The Irwin-Hall distribution (sum of standard uniform) offers this characteristic, but is not treating non-integer orders $N$. However, if you e.g. simply independently sum up e.g. 2 standard uniform $[0,1]$ and one 3rd with a smaller range like $[0,0.25]$ you will indeed nicely obtain a more general and smoothly extended version of Irwin-Hall for any arbitrary order (like $N=2.25$ in this case). However, I wonder if it is possible to find a practical closed formula for the CDF?


title: T-statistic of correlation coefficient

text: Hi I am trying to calculate the t-statistic for a correlation coefficient between two vectors $x$ and $y$. The individual vectors shows signs of autocorrelation. I have made use of the formula: $\frac{r}{\sqrt{(1-r^2)(n-2)}}$ but I'm not sure this is the correct way of solving the problem. The correlation between the two vectors has been calculated with the Excel formula for correlation.


title: Analyzing Difference in Change between two groups with unpaired samples

text: I'm currently working with a data set where there is a control group and an intervention group. There are a couple of data points where we're interested in seeing if the intervention group has a statistically significant difference in effect (basically is $x_2 - x_1$ in intervention different than $x_2-x_1$ in control). Unfortunately, the samples collected within the groups are not paired - in other words, the same people were not interviewed at baseline and endline. I can use an unpaired t-test to compare endline-to-endline, or baseline-to- endline for the groups individually. Were the samples paired, I could just take the difference and use an unpaired t-test to compare the differences. Since they're unpaired, I'm at a bit of a loss. Is there a way to adjust the t-test for this? If not, are there appropriate statistical tests that can be done for this situation without artificially pairing up the samples? (As is too commonly the case, I was not involved in the design of this experiment, but I hope I've conveyed to the group the importance of paired sampling for the future.)


title: Modelling relation between two persistent AR(1) processes? Is my approach reasonable and how to adjust standard errors?

text: I want to run a simulation in which I want to find out whether there is a relation between the independent variable $x_t$ and the dependent variable $y_t$. I.e., in the following regression I want to find out if $\beta$ is signficantly different from zero: $$ y_t = \beta x_t + \epsilon_t. $$ In empirical data, $y_t$ is stationary, but close to a random walk.I want to model it as an AR(1) process. Of course, I could just model $x_t$ as an AR(1) process and if $\beta$ is different from zero, $y_t$ would end up as such a process as well. However, I also want to consider the case in which $\beta=0$ and I don't want to confuse the reader with changing descriptions for different scenarios (i.e., "if $\beta=1$ I simulate $y_t$ from the above equation, if $\beta=0$, I simulate $y_t$ as an AR(1) process.") I came up with the following solution to this dilemma. Why not model both $y_t$ and $x_t$ as AR(1) processes with potentially correlated errors. That is: $$ y_{t+1} = \tau_y y_t + u_{t+1}, \quad 0 < \tau_y < 1, \quad u_{t+1} \sim N(0, \sigma_u^2) $$ and $$ x_{t+1} = \tau_x x_t + v_{t+1} , \quad 0 < \tau_x < 1, \quad v_{t+1} \sim N(0, \sigma_{v}^2), $$ where $$ cov\left(\begin{bmatrix} v_{t+1} \\ u_{t+1} \end{bmatrix}, \begin{bmatrix} v_{t+1} & u_{t+1} \end{bmatrix}\right) = \begin{bmatrix} \sigma_v^2 & \sigma_{v} \sigma_u \rho_{u,v} \\ \sigma_{v} \sigma_u \rho_{u,v} & \sigma_u^2 \end{bmatrix} $$ Since $\beta = \frac{Cov(y, x)}{\sigma_x^2}$ and $$ \begin{aligned} Cov(y, x) &= E[y_{t+1}x_{t+1}] \\ &= E[(\tau_y y_t + u_{t+1})(\tau_x x_t + v_{t+1})] \\ &= \tau_y \tau_x E[y_t x_t] + E[u_{t+1} v_{t+1}] \\ &= \tau_y \tau_x Cov(y, x) + \sigma_{u} \sigma_{v} \rho_{v, u} \\ \end{aligned} $$ we get by rearranging $$ \beta = \frac{\sigma_u}{\sigma_{v}}\rho_{v, u} \frac{1 - \tau_x^2}{1-\tau_x \tau_y}. $$ Fair enough, this is a little bit complicated, but it allows me to control $\beta$ without changing the structure of either $y_t$ or $x_t$ and both those processes are AR(1) in this setup. I basically have two questions:   1. Is this a reasonable approach given the requirement that preferably both $x_t$ and $y_t$ should be AR(1) processes. The interpretation of the relation should really be just like in a simple OLS setup, i.e. has $x_t$ an impact on $y_t$ (I'm not concerned about causality here, just relation). Would you consider my setup still as a reasonable way of modelling it? (I just don't make a direct link, but use the error structure for that. I don't see a problem with this approach. If $x_t$ and $y_t$ are related via correlated shocks, so be it.)   2. I already simulated this and it works fine. However, I also want to obtain correct confidence intervals for each run. That is, I want to run the regression in the first equation, get $\widehat{\beta}$ and also standard errors that should be reasonable. Now I noticed that the errors are highly autocorrelated and this autocorrelation seems to be identical to $\tau_y$. However, I could not formalize it. So it would be great to know how I would have to adjust the standard errors. (It would be extra awesome if I could get a hint how to implement that in R.)              set.seed(123)     library(MASS)     ### Set start values     nrT <- 1e5     burnin <- 1e3     sd_y_shock  <- 1     sd_x_shock  <- 1     corr_x_y    <- 0.5     tau_y  <- 0.8     tau_x  <- 0.1     ### Simulate the correlated shocks     shocks <- mvrnorm(nrT + burnin,                        mu = c(0,0),                        Sigma = matrix(c(sd_y_shock^2, corr_x_y * sd_y_shock * sd_x_shock,                                         corr_x_y * sd_y_shock * sd_x_shock, sd_x_shock^2),                                       nrow=2))          vec_y <- arima.sim(list(order = c(1, 0, 0), ar = tau_y),                         n = nrT + burnin,                          innov = shocks[, 1])           vec_x <- arima.sim(list(order = c(1, 0, 0), ar = tau_x),                        n = nrT + burnin,                         innov = shocks[,2])     ### Check that formula derived above is correct; the two should be similar     sd_y_shock/sd_x_shock * corr_x_y * (1 - tau_x^2)/(1 - tau_x * tau_y)      coef(lm(vec_y ~ vec_x))[2]     ### Plot ACF     # Note that, independent from corr_x_y and tau_x, the autocorrelation structure seems     # always to be tau_y     acf(lm(vec_y ~ vec_x)$resid)      ![enter image description here](http://i.stack.imgur.com/74Jjm.png)


title: Need to find cumulative probability distribution given this data set

text: So I have a data set like so:               Email                  Posted           Received                            Platform     ss104@gmail.com 3/20/2012 1:11  Tue, 20 Mar 2012 01:07:05 -0700 (PDT)   Wikipedia     ss104@gmail.com 3/20/2012 1:11  Tue, 20 Mar 2012 01:07:22 -0700 (PDT)   Mailing list     ss104@gmail.com 3/20/2012 1:11  Tue, 20 Mar 2012 09:08:27 -0700 (PDT)   Twitter     ss104@gmail.com 3/20/2012 1:11  Tue, 20 Mar 2012 12:22:46 -0700 (PDT)   Wikipedia     ss104@gmail.com 3/20/2012 1:11  Tue, 20 Mar 2012 17:45:51 -0700 (PDT)   Mailing list     ...      What this means is that a certain email address was posted at a certain time on a certain platform. It received an email at the time under `Received` column. What I want to find is which platform sends the most amount of emails per week. (I'll be using python and matplotlib for visualizing). CPD (cumulative probability distribution) for this data is a good start. I'm not sure how to start with this. Do I just count the occurrence of each platform, and divide it by the summation of difference in times (`Received - Posted` in weeks) for each email of that platform? And then I could also make histogram or box plot for the data, but I'm not what actually to compute. Can anyone give me pointers?


title: How to show that polar coordinates in a uniform distribution on a disk are independent?

text: Let the random point $(X,Y)$ be uniformly distributed on the unit disc $D=\{(x,y):x^{2}+y^{2}<1\}$. Show that the polar coordinates $R\in [0,1)$ and $\theta \in [0,2\pi)$ of the point are independent. Can you help me with this exercise please?


title: How to describe variation on multiple levels?

text: This is continuation for a series of questions (1, 2). I have a data set from an experiment with 2x2 design. A replicate consists of 20-50 normally distributed replicate measurements. Each treatment combination has 3 of these replicates. Measured average response for each replicate seems not to be normally distributed. ![enter image description here](http://i.stack.imgur.com/tfJEe.png) Figure 1. Example of the problem. Error bars show standard error of mean. Since this is an experiment and I am supposed to talk about an effect of a treatment, I would like to combine the variation within measurements (error bars) and between replicates (dots in the figure) into one number / treatment combination and use the error bars with one mean value to show the overall effect of a treatment (a dot instead of three). Ideally I would like to use confidence intervals of some sort (asked here). I have understood that since this is an unbalanced nested design, I should use mixed models or a simple 2-way ANOVA for mean values, but how to compress the variation on two levels into one figure? Would the delta method be something for me? R code with an example of the data:               library(ggplot2)          x1 <- c(rnorm(50, 14,3),rnorm(35, 7,1),rnorm(40, 15,9))     x2 <- c(rnorm(43, 6,3),rnorm(32, 7,1),rnorm(40, 8,4))     x3 <- c(rnorm(50, 15,5), rnorm(50, 10,7), rnorm(50, 13,9))     x4 <- c(rnorm(26, 14,2), rnorm(43, 25,10), rnorm(45, 15,9))          dx1 <- data.frame(Treatment = rep("T1", length(x1)), Temp = rep(10, length(x1)), Rep = rep(c(1,2,3), times = c(50,35,40)), Meas = x1)     dx2 <- data.frame(Treatment = rep("T1", length(x2)), Temp = rep(20, length(x2)), Rep = rep(c(4,5,6), times = c(43,32,40)), Meas = x2)     dx3 <- data.frame(Treatment = rep("T2", length(x3)), Temp = rep(10, length(x3)), Rep = rep(c(7,8,9), times = c(50,50,50)), Meas = x3)     dx4 <- data.frame(Treatment = rep("T2", length(x4)), Temp = rep(20, length(x4)), Rep = rep(c(10,11,12), times = c(26,43,45)), Meas = x4)          # Entire data set          dat <- rbind(dx1,dx2,dx3,dx4)          # Plot overview          w <- ggplot(dat, aes(x = factor(Rep), y = Meas))     w + geom_boxplot(aes(fill = factor(Temp)))          # Averages (with se)          aveg <- aggregate(Meas ~ Treatment + Temp + Rep, data = dat, FUN = mean)     se <- function(x) sd(x)/sqrt(length(x))     SE <- aggregate(Meas ~ Treatment + Temp + Rep, data = dat, FUN = se)          dat2 <- merge(aveg, SE, by = c("Treatment", "Temp", "Rep"), sort = F)          colnames(dat2)[colnames(dat2) %in% grep("\\.x", colnames(dat2), value = T)] <- "mean"     colnames(dat2)[colnames(dat2) %in% grep("\\.y", colnames(dat2), value = T)] <- "se"          # Plot entire data set          p <- ggplot(dat2, aes(x = Treatment, y = mean, ymax = mean + se/2,      ymin = mean - se/2))          p + geom_pointrange(aes(color = factor(Temp)),      position=position_dodge(width=0.50), size = 1)


title: How actually my time series models tells me about data behaviour?

text: I have different models namely AR, ARMA, ARIMA for different time series (in data set of 4449) stationary-3096, Non stationary-1353 , why this classification appears to be? If I have $AR(p=1-5)$ and $AR(p >5)$ what does my data tell ?, I have many time series which follows ARMA-1239 best models, some are AR-90 best model & ARIMA-7 best models out of 3096 stationary data set Please, kindly tell me in philosophical/non-statistical way (because I didn't understand by term like normal/marginal distribution).


title: Algorithm for finding all local maxima of a boosted regression tree?

text: Given a weighted sum of regression trees, is there an efficient algorithm to find **all** local maxima ? (I would tend to think that a gradient based method will find **some** maxima but it is not entirely clear to me how to ensure that one will find all of them)


title: When are Bayesian methods preferable to Frequentist?

text: I really want to learn about Bayesian techniques, so I have been trying to teach myself a bit. However, I am having a hard time seeing when using Bayesian techniques ever confer an advantage over Frequentist methods. For example: I've seen in the literature a bit about how some use informative priors whereas others use non-informative prior. But if you're using a non- informative prior (which seems really common?) and you find that the posterior distribution is, say, a beta distribution...couldn't you have just fit a beta distribution in the beginning and called it good? I don't see how constructing a prior distribution that tells you nothing...can, well, really tell you anything? It does turn out that some methods I have been using in R use a mixture of Bayesian and Frequentist methods (the authors acknowledge this is somewhat inconsistent) and I cannot even discern what parts are Bayesian. Aside from distribution fitting, I can't even figure out HOW you would use Bayesian methods. Is there "Bayesian regression"? What would that look like? All I can imagine is guessing at the underlying distribution over and over again while the Frequentist thinks about the data some, eyeballs it, sees a Poisson distribution and runs a GLM. (This isn't a criticism...I really just don't understand!) So..maybe some elementary examples would help? And if you know of some practical references for real beginners like myself, that would be really helpful too!


title: Testing a model with truncated data

text: This is probably a simple question. I'm studying events which have N outcomes, of which exactly one is correct. N is very large, more than a billion (and is known). There are many possible events, some of which are tested multiple times. I would like to test the following model: a given event is tested correctly with probability 1 - p. Otherwise, the test result is chosen uniformly at random from the N outcomes (hence, almost certainly incorrect). I want to (1) see if this model is plausible and (2) if so, find p. p should be small, just a few percent. The data I have omits events where two or more tests have the same value. (Because N is large, if two tests have the same value that is likely to be the correct result.) So what I have is, for each k, the number of different events which have had k different tests, no two of which match. Most of the events have had k = 1 tests but some have 2, 3, or even more tests. In most cases one result is correct and the rest are incorrect; in some cases all are incorrect. How can I test my hypothesis?


title: Software for LASSO for high dimensional dataset

text: I have a dataset with between 10,000 and 100,000 feature values. The number of datapoints is between 1,000 and 10,000. I want to perform a LASSO on this dataset but can't really find any good software to do so. Does anyone have any suggestions?


title: How to write a function to generate a sequence of points in R?

text: This is the PDF that I am dealing with:               fx = 0.3 if (0<=x<1)      0.1 if (1<=x<2)      0.25 if (2<=x<3)     0.15 if (3<=x<4)     0.2 if (4<=x<5)      0 otherwise      I have to write a function `gen_xy` that will generate a sequence of points `(X,Y)` uniformly distributed in `(0,5)x(0,0.5)` until one lies in the region under the curve of `y=fx(x)` This is my code so far:                   fx <- function(x) c(0, 0.3,0.1,0.25,0.15,0.20, 0)         [findInterval(x, c(-Inf, 0:5,  Inf))]              x <- runif(n,0,5)         fx <- stepfun(x = 0:5, y = c(0,0.3,0.1,0.25,0.15,0.20,0))         plot(fx, ylim = c(0,0.5),xlim = c(0,5), main = 'f(x)')      Now this is my attempt at writing the code for my function:                   gen_xy <- function()         {          done=0          while(done==0)          {           x=runif(1,0,5)           y=runif(1,0,0.5)           print(c("x,y",c(x,y)))           if(y < fx(x))           {            done=1           }          }          xy=c(x,y)          xy         }      But I think the part `if(y < fx)` is wrong? I then need to generate a sample of 1000 points and plot them to check if they are an appropriate sample from under the curve of `y=fx`. How would I go about writing such a code?


title: What is covariance matrix adaptation evolution strategy?

text: I am new to evolutionary algorithm. I have studied Covariance Matrix Adaptation Evolution Strategy. I am not good at statistics. So could you please explain me in simple language (I mean not too many equations)   1. What is CMA-ES?   2. How does it work?   3. Why is it superior to other strategies?


title: Multiclass logistic regression with mlogit in R

text: I have a multi-class dataset like the following (a,b,c,d are features and e is the class (it can be 0,1 and 2)).                   a b c d e     1   1 1 2 2 1     2   1 2 4 2 0     3   1 2 4 2 0     4   2 2 2 2 0     5   2 1 2 2 2      I am trying to use mlogit package in order to see which column is more important but I am having a difficulty to understand how to use it. What I am doing is:               dataset$e<-as.factor(dataset$e)     mldata<-mlogit.data(dataset, choice="e")          > mldata[1:5,]         a b c d     e chid alt     1.0 1 1 2 2 FALSE    1   0     1.1 1 1 2 2  TRUE    1   1     1.2 1 1 2 2 FALSE    1   2     2.0 1 2 4 2  TRUE    2   0     2.1 1 2 4 2 FALSE    2   1      Now, in order to see the coefficients, I am constructing the model like that: > mlogit.model<- mlogit(e~1|a+b+c+d, data = mldata, reflevel="1") mlogit.model               Call:     mlogit(formula = e ~ 1 | a + b + c + d, data = mldata, reflevel = "1",     method = "nr", print.level = 0)          Coefficients:          alt0       alt2     alt0:a     alt2:a     alt0:b     alt2:b     alt0:c       -211.0953   -89.7558    27.9911    33.1440    37.7503     7.3585     8.5072          alt2:c     alt0:d     alt2:d          3.0950    37.1340     5.5584        But now I don't understand what are alt0, alt2? What are the real coefficients of a,b,c and d?


title: Identifiability and unbiasedness

text: How do you show that if my model parameter $\theta$ (scalar) is U-estimable (i.e. if there exists an unbiased estimator of $\theta)$, then $\theta$ is identifiable? This makes sense intuitively, but I'm not sure how to prove it.


title: Inverse prediction of percentile

text: Suppose we consider a linear normal regression $y\sim$Normal($\mu,\sigma$) and $\mu=a+bx$. I have seen documentations for methods for obtaining a confidence interval for $x$ for a specific mean value of $y$. Suppose instead of the mean $y$, I want a confidence interval in $x$ for top 5%ile of $y$ (e.g., so that prediction in $x$ will be lower than that for the mean $y$, provided $b>0$). I found a method for getting a confidence interval of percentile in $y$ (https://onlinecourses.science.psu.edu/stat414/node/317). But I don't know what a next step is. Ideally, I want a method that works for any other distributions (e.g., gamma distribution and log-normal distribution).


title: How do I get better forecasts for this data

text: How do I get better forecasts for my model? Is the plot supposed to look like this? I am using the code:               fit <- auto.arima(blah)     fcast <- forecast(fit,100)      ![enter image description here](http://i.stack.imgur.com/TeiGU.png)


title: How to model repeated measures' results in meta-analysis?

text: I am developing network meta-analysis BUGS model that I will use to compare Hazard Ratios in studies that assume Box Proportional Hazard model. As input I am given pairs: (log HR, its standard error) that measure magnitude of the effect of difference between two treatments. Naturally I assume normal distribution model for the log HR ("LHR" for short) with its given SE. Problem is, when I have more than one LHR result from a single study: sometimes these results are based on independent groups - and then I treat them just as if they came from separate studies. But more often than not these results are based on multiple comparisons in more than two-arms experiment on the same group. It means, that they have e.g. 4 groups of patients: A, B, C and d and report Hazard Ratios A vs B, A vs C, A vd D, B vs C, B vs D and C vs D. These 6 results are not-so-independent, so I can't really treat them as independent in my BUGS model. So I believe that treating them as separate studies is wrong, but _what else can I do_? If that is relevant for you, I can post actual BUGS model.


title: Is a pre-post intervention study with a control group without random assignment quasi-experimental?

text: I am currently rather confused about the design type for my study. I am looking at a number of measures pre and post intervention with control group in a community sample. The sample collected was not through randomization. Is the design type quasi-experimental? What kind of statistics can i use for this design. Note: Sample size is 20 in each group before and after.


title: Logistic Regression – Evaluation Metrics

text: I have more of a programming background, and I am fairly new to statistics. I am currently trying to solve some sample exercises to get more familiar with data science / modelling. ### Problem Background A user posts a request on a forum. Considering number of responses / number of up votes / users' reputation / etc., can you predict if the request will be fulfilled? There are about 4000 rows of test data available. ### Current Attempts I have created basic graphs to investigate correlation, and applied a few transformations to the data. I have decided to go with logistic regression, and based on the graphs, have chosen 2 factors for the initial model. ### Problem I have created the basic model using R, and I can see some summary stats, but I am a bit lost when I try to understand (a) How well the model fits the data and (b) whether the model is accurate. I tried Googling around, but most of the articles were too technical for me. Is there a cheat sheet/quick test of some sort that I can apply? Or can somebody suggest the simplest "complicated" article?


title: How to read a value-quantile diagram?

text: I've been reading an article on ACM queue, and got confused by figure 6 (latency distribution of a push-only workload), where the authors used the y axis to display the latency, and the x axis to denote the quantile. ![enter image description here](http://i.stack.imgur.com/ptBuo.png) The problem is, I really don't get how to read that plot:   * What do the bars show? How can interpret a long bar at 0.95, for instance?   * Where do I find, e.g. the mean? By the way, what is the name of this type of plots?


title: Hierarchical text classification in R

text: I am working Automatic Email routing. My historical data set has Email Description,Division and Category. There are around 20 division and 150 categories. I have built one vs. all classifier (RTextTools - SVM) for 20 divisions and trained 20 classifier and when the new email comes in , my container(SVM Models) gives me predicted Division (with good accuracy). My question is from here( after getting Division, since i have got good accuracy level using one vs. all approach - ) - how would i predict Categories . As there are 150 levels - i cannot do one vs. all - as did for Division.And more than 80% of data the is covered with in 40 Category levels. Can somebody suggest a way of getting Category level prediction(second level) - or any reference how it can be done in R.


title: Incorporating function of random variable in bayesian network.

text: Let A be a random variable hence a node in my bayesian network. I wish to apply some function on A and then connect it to some other node in the network. So how to do that?


title: Picking socks probability proof

text: I'm trying to follow the proof from question 1b out of the book "50 Challenging Problems in Probability" by Mosteller. The problem states: A drawer contains red and black socks. When two socks are drawn at random, the probability that both are red is $1/2$. How small can the number of socks in the drawer be if the number of black socks is even? The answer given starts as follows: Let there be $r$ red and $b$ black socks. The probability of the first socks' being red is $\frac{r}{r+b}$ and if the first sock is red, the probability of the second's being red now that a red has been removed is $\frac{r-1}{r+b-1}$. Then we require the probability that both are red to be $\frac{1}{2}$, or $$\frac{r}{r+b} \times \frac{r-1}{r+b-1} = \frac{1}{2}$$ Notice that $$\frac{r}{r+b} > \frac{r-1}{r+b-1} \quad \text{for b} >0$$ Therefore we can create the inequalities $$\left(\frac{r}{r+b}\right)^2 > \frac{1}{2} > \left(\frac{r-1}{r+b-1}\right)^2$$ This is where I'm confused. Why is it that $\left(\frac{r}{r+b}\right)^2 > \frac{1}{2}$? If $r=1, b=100$ then obviously $\left(\frac{1}{101}\right)^2 < \frac{1}{2}$. Am I missing some obvious assumption?


title: Is a restricted Boltzmann machine trained using a contrastive divergence algorithm an unsupervised learning method?

text: I am designing RBM with CD training methodology in MATLAB. CD is an unsupervised learning algorithm, so I do not need to use any label information included in input. Does that make sense?


title: In a linear regression, how do I find the predictor for which the response is zero?

text: I'm linearly regressing some response `y` onto some predictor `x`. I'm interested in knowing for what `x` does `y = 0`. I can think of two ways to do this. Let me illustrate with some sample data:               x <- 1:10     y <- 20 - 2 * x + rnorm(10)      I can either linearly regress `y` onto `x` and solve the equation explicitly:               - coef(lm(y ~ x))[1] / coef(lm(y ~ x))[2]     (Intercept)         10.29915       Or I can try to be clever and observe that my problem is equivalent to regressing `x` onto `y` and predicting `x` for `y = 0`:               coef(lm(x ~ y))[1]     (Intercept)         10.19658      However, this doesn't give quite the same result. So which approach is correct?


title: Wilcoxon test in boot() function

text: I'm trying out the boot() function for internal validation of a logistic glm model using the AUC (aka c-statistic) as my performance measure. My problem is that depending on the dataset I use, sometimes the function gets "50 or more warnings," all of which are of this type:               Warning messages:     1: In wilcox.test.default(pred[obs == 1], pred[obs == 0],  ... :       cannot compute exact p-value with ties      Even though the function produces these warnings, it does also produce the relevant bootstrap statistics. I'm not sure what to make of this situation since I don't even know in what context boot() is calling wilcox.test. The documentation ?boot doesn't mention the use of the Wilcoxon test. Here is an example that generates these warnings:               data(nuclear)     library(verification)          AUC = function(data, i) {         d = data[i,]         rs1 = glm(ne ~ cost, family=binomial(link="logit"), data=d)         return(roc.area(d$ne, predict(rs1))$A)     }          library(boot)     boot(data=nuclear, statistic=AUC, R=600)      And an example that doesn't:               data(melanoma)          AUC = function(data, i) {         d = data[i,]         rs1 = glm(ulcer ~ thickness, family=binomial(link="logit"), data=d)         return(roc.area(d$ulcer, predict(rs1))$A)     }          boot(data=melanoma, statistic=AUC, R=600)      (I apologize if this question belongs in Stack Overflow instead of here -- not sure how statistical vs. programmatic the issue is.)


title: When is it appropriate to use a paired Wilcoxon test?

text: Could someone please explain in plain English (preferably without using the word paired) or by using examples when it is appropriate to use `wilcox.test(x,y,paired=TRUE)`? I have read the R documentation and various links on google and I'm not 100% sure that I get it. If someone would be kind enough just to explain it simply I would feel much happier going ahead and using it on my data. At the minute I think it's appropriate to compare data like `amount.on.monday` and `amount.on.friday` but not for comparing `time.i.eat` and `amount.on.monday`?               time.i.eat=c(1,2,3,4)     amount.on.monday=c(1,1,1,1)     amount.on.friday=c(2,2,2,2)


title: Creating power law distributed series in R

text: I would like to create a power law distributed series with 100 numbers and a minimum value of 10. I could not find any command in R that directly does this. What is the easiest way?


title: Some clear range and interval definitions

text: What is your or an offical (please provide link to citation) definition of range AND interval **___**? Or perhaps put in another way: what is the major important difference between the two terms **___**? My perception is the following: > Statistical definition of "range" is between some defined data points - like > max and min values (allows you to state that all data is within the range of > max and min value) AND > Statistical definition of interval is used to put equal distance(s) above > and below OR between specific data OR calculated value(s). So example a 95% > confidence interval around the mean means that 95% of all observation are > expected to occur within this interval. Another example would be a 10 % > interval between the min and max which will create 9 equally distanced > interval points between the min and max value. Do you have other definitions AND/OR examples of ranges AND/OR intervals **___ __ ___**?


title: Singular covariance matrix in exploratory Factor Analysis

text: I'm kind of a noob to EFA and am trying to use the FANode object in Python. This is from the MDP library. I am using it on survey data to see which variables are tied together. Whenever I run it on my data, I get the following error:               mdp.NodeException: The covariance matrix of the data is singular. Redundant dimensions need to be removed.      I was wondering if anyone else has experienced these errors? What do they mean to you? Again I am quite the noob and would prefer to ask the community before diving into a text book.


title: lmer giving worse performance than lm

text: I am having trouble training a model for nested data about house prices. Lets say my data looks like following:                 logPrice bedCount bathCount                city      0.6517920        4       2-3        Redwood City      0.4402192        1       1-2 South San Francisco      0.5922396        2       1-2           San Mateo      0.4606918        3       1-2 South San Francisco      0.7592523    5plus       3-4           San Mateo      0.4710397        1       1-2        Redwood City      `bedCount`, `bathCount` and `city` are factors. As a baseline, I trained a simple linear model ignoreing nested structure of the data (houses are nested within cities).               lm.model = lm(formula = logPrice ~ 1 + bedCount + bathCount + city)      which corresponds to following assumption: `logPrice`$_i = \beta_0 + \beta_1\cdot$ `bedCount`$_i + \beta_2\cdot$ `bathCount`$_i + \beta_{3,j[i]}\cdot I$(`city`$_{j[i]}) + \epsilon_i$ where $\epsilon_i \sim N(0, \sigma^2_{logPrice})$ and $I$(`city`$_{j[i]}$) is the indicator function for city of the $i^{th}$ house (which is 1). Now, I trained a 2-level hierarchical model:               lmer.model = lmer(formula = logPrice ~ 1 + bedCount + bathCount + (1 | city))      which corresponds to the following assumption: `logPrice`$_i = \beta_0 + \beta_1\cdot$ `bedCount`$_i + \beta_2\cdot$ `bathCount`$_i + \beta_{3,j[i]}\cdot I$(`city`$_{j[i]}) + \epsilon_i$ where $\epsilon_i \sim N(0, \sigma^2_{logPrice})$ and $\beta_{3,j} \sim N(0, \sigma^2_{\beta_3})$ Now, on the training data, `lm.model` gives me lesser average RMSE than `lmer.model` which shouldn't happen because linear regression is a special case of multilevel linear regression (I didn't care to check average RMSE on test data because that on training data itself should be lower for 2nd model than that for 1st model). In fact, my data has multiple levels (houses nested within subdivisions, which are nested within zipcodes, which in turn are nested within cities) and the performance gets worse as I add more and more levels to the model (i.e. model with random effect `(1 | subdivision)` does worse than that with random effect `(1 | zipcode) + (1 | zipcode:subdivision)`, which in turn does worse than a model with random effect `(1 | city) + (1 | city:zipcode) + (1 | city:zipcode:subdivision)`). What am I missing?


title: Computing the F-ratio under null hypothesis

text: Statistical model for a Completely Randomized Design: $$y_{ij}=\mu+\tau_{i}+\epsilon_{ij}\quad i=1,2,\ldots,a; j=1,2,\ldots,n$$ where, $y_{ij}$ is the $j$th observation of $i$th treatment effect $\mu$ is the overall mean $\tau_{i}$ is the $i$th treatment effect $\epsilon_{ij}$ random error term and $\epsilon_{ij}\sim NID(0,\sigma^2)$ For testing the hypothesis of no difference in treatment means $$H_0:\tau_1=\tau_2=\ldots=\tau_a = 0 \quad vs.\quad H_1:\tau_i\ne 0\quad \text{for some } i.$$ Under `null` hypothesis , the ratio $$F=\frac{\frac{SS_{treatments}}{a-1}}{\frac{SS_{error}}{N-a}}$$ is distributed as F with $a-1$ and $N-a$ degrees of freedom .   * My question is: Why do i have to compute the F-ratio under `null` hypothesis, why do not under `alternative` hypothesis ?


title: Actual value falling within prediction interval rather than predicted value falling within actual value range?

text: I'm trying to predict a price for an item on the website. I'm interested in classifying the prediction as good or bad. There are two approaches that I can think of, and I'm wondering which one is better and why:   1. See what percent of the time the predicted value falls within a certain distance away from the actual value. If it falls within a certain interval of the actual value, classify this as good.   2. See what percent of the time the actual values fall within a prediction interval of the predicted value. If the actual value falls within the prediction interval, classify as good. Let $in$ be the event that the actual value falls within the prediction interval. The prediction interval is defined as $\hat{\mu} \pm x*\frac{\sigma}{\sqrt{n}}$ where $\hat{\mu}$ is the predicted value, $\sigma$ is the variance of the sample, and n is the number of observations. $x$ is a value chosen to represent the level of desired certainty. Let A represent the actual value: $P(in)=P(\hat{\mu} - x*\frac{\sigma}{\sqrt{n}} < A < \hat{\mu} + x*\frac{\sigma}{\sqrt{n}})$ Which is the correct approach, and why? * * * *Note: I have since realized that I may have incorrectly termed the interval displayed above a prediction interval when in fact it is a confidence interval. What I would like to use is the prediction interval...still reading up to understand the difference between these intervals.


title: Is it possible to accept the alternative hypothesis?

text: I'm aware of several related questions here (e.g., Hypothesis testing terminology surrounding null, Is it possible to prove a null hypothesis?) but I don't know the definitive answer for my question below. Suppose a hypothesis test where we want to test whether a coin is fair or not. We have two hypotheses: $H_0: p(head)=0.5$ $H_1: p(head)\neq0.5$ Suppose we use 5% significance level, there are two possible cases:   1. When we obtain the data and find that the p-value is less than 0.05, we say "With significance level 5%, we reject $H_0$."    2. p-value is greater than 0.05, then we say "With 5% significance level, we cannot reject $H_0$." My question is: > In case 1, is it correct to say "we accept $H_1$"? Intuitively, and from what I have learned in the past, I feel that "accepting" anything as a result of hypothesis testing is always incorrect. On the other hand, in this case, since the union on $H_0$ of $H_1$ covers the whole "space", "reject $H_0$" and "accepting $H_1$" look exactly the same to me. On another thought, I can also think of the following idea, which says it is incorrect to say "we accept $H_1$": > We have an evidence strong enough to believe that $H_0$ is not true, but we > may not have an evidence strong enough to believe that $H_1$ is true. > Therefore, "rejecting $H_0$" does not automatically imply "accepting $H_1$" So, what is the right answer?


title: Probability of $k$ successes in no more than $n$ Bernoulli trials

text: I feel like these questions should be related to the binomial and geometric distributions, but I can't quite work out how. Suppose I do sequential independent Bernoulli trials each with success probability $p$. Rather than fixing the number of the trials, as with the binomial distribution, I want to fix the number of successes I want to achieve. Call this number $k$. Once I obtain $k$ successes, I stop the Bernoulli trials. What is the probability that I perform exactly $n$ trials? What is the probability that I need to do _no more_ than $n$ trials to achieve $k$ successes? As a concrete example, suppose I'm flipping a fair coin and I will continue to do so until I see 20 heads. What is the probability that I see the 20th head on the 40th coin toss? What is the probability that I have to do 40 or fewer tosses until I see 20 heads? Edit: Ha, I don't know what it is about the answer coming to you as soon as you give up and ask for help, but I've got it now. The probability of needing exactly $n$ trials for $k$ successes is just the binomial probability of $k−1$ successes in $n−1$ trials times the probability of success on the $n$th trial. The probability of needing no more than $n$ trials is the sum of the probabilities of needing exactly $k$ through $n$ trials. I'm pretty sure that works, let me know if you think this is wrong. That said, if anyone has any ideas for turning this into a nicer, more compact expression, I'm all ears!


title: Is this the definition of over-fitting?

text: Overfitting is when we have a model which has memorized the training data and does not perform well in real-world cases. Okay, say that I had some training points which look like this: ![enter image description here](http://i.stack.imgur.com/NEzOx.png) What if the red curve was the actual 'real-world' relationship. And I found this exact model through a learning algorithm based on those training observations. Has my model overfit by definition despite it being the 'real-world' relationship? I am assuming yes, but I just want to make sure. Thanks


title: Data normalization before giving to Neural Nets or Deep Learning algorithm?

text: What kind of normalization scheman is required for the best of NN algorithms? I saw some people just give the data to signum function before passing to NN and some of those process data by regular standardization. What is the best result method?


title: category selection with LASSO

text: Suppose one has two features: color = {R, G, B} and t-shirt size = {S, M, L} and wants to regress these features on the probability of a sale, call it p. So the model is p ~ color + size. Now, the design matrix will have an intercept, and colors can be represented as               R -> [1 0 0]     G -> [1 1 0]     B -> [1 0 1]      Similarly for sizes. But color _and_ size together will have 3x3=9 combinations, represented as (color first, then size):               R + S -> [1 0 0 0 0]     R + M -> [1 0 0 1 0]     R + L -> [1 0 0 0 1]     G + S -> [1 1 0 0 0]     G + M -> [1 1 0 1 0]     G + L -> [1 1 0 0 1]     B + S -> [1 0 1 0 0]     B + M -> [1 0 1 1 0]     B + L -> [1 0 1 0 1]      Now, when running LASSO on this, R + S are mixed together, i.e. say all coefficients but the intercept get shrunk to zero, which means that the likelihood is dependent on the color R and t-shirt size S, but how can one tell if R and/or S is more important? Is there another way to represent these categorical predictors without mixing two categories from different features in the intercept? Or is one stuck repeating the feature selection on R + S in this case to tell which one is more important? Is LASSO even the right algorithm to use for this setting, or would one just resort to classical backward/forward feature selection and then category selection using BIC or cross-validation as the error metric?


title: Mean of means -> error propagation or uncertainty or both?

text: my problem is as follows. I have a simulation of a neural network which creates activity patterns, learns them and then tries to retrieve previously learned patterns one by one. The performance of these retrievals (=remembering) is the result I'm interested in. Simplified, this is how it works (I used fictive numbers to make it more understandable):   * Create 40 random activity patterns within a network of 200 neurons.   * Learn those patterns.   * Retrieve (=remember) each pattern. This gives me a retrieval performance (between 0 = "not at all" and 1 = "perfectly"). Now, this retrieval process uses randomly generated numbers and will be performed 40 times for each pattern leading to a mean performance for the given pattern set. I also calculate sample standard deviation and uncertainty of the mean. My problem is now this: I change some parameters, say the number of patterns created in the first step, and for each number of patterns, I run the whole thing (create patterns, learn them, retrieve them) 40 times. Then, I have for each number of patterns 40 mean performance values plus their corresponding standard deviations and uncertainties of the mean. I know want a single mean (of all these 40 performances) for each number of patterns. That's easy, but how do I calculate the uncertainty of this overall performance value? Should I use gaussian error propagation for the uncertainties of the mean values I'm averaging over? Or should I just use the uncertainty of the newly calculated overall mean performance?? So, this is how it looks like:               # of patterns    (mean performance)1    (mean performance)2    ...     10               0.82 +- 0.01            0.79 +- 0.05          ...     11               0.75 +- 0.02            0.77 +- 0.03          ...     12               ...                     ...     ...      I calculate the mean for all (mean performance) values corresponding to a number of patterns, so basically I'm taking the mean of mean values. If you could give me a hint of how to best calculate the overall uncertainty of the result, I'd very much appreciate it. I've googled a lot and also read some information sheets about errors in measurements but am not sure as how to best proceed. My current idea is to just use gaussian error propagation but then I won't take the uncertainty of the newly calculated mean value into account. Thanks for your help! Daniel


title: Sampling a high dimensional function

text: I have a function $f$ that takes $N$ real-valued inputs and is very expensive to compute. I am studying the properties of $f$ by generating samples. Each input $n$ has a range of values $(n_{min}, n_{max})$ that it can take on. Let's say that due to time constraints I can only generate 1000 samples. Is it "better" to select the set of $N$s that I feed to $f$ by (A) iterating from $n_{min}$ to $n_{max}$ with a big enough step size for each $n$, or (B) uniformly sampling each $n$ on the range of $(n_{min}, n_{max})$? Choice (A) has the desirable property of keeping every other input fixed while varying only one value at a time, but choice (B) has the desirable property of likely exploring more parts of the input space.


title: Multiple Bernoulli and Multinomial Distirbution

text: It's well known that language can be modeled by Multinomial distribution and Multiple Bernoulli distribution. So far I don't see any advantage of Multiple Bernoulli distribution representation over Multinomial representation. Both models equal computationally, Multinomial considers the number of occurrence of the unigram/bigram/.. tehrefore it's more precise. What's the case when Multiple Bernoulli is preferable over Multinomial.


title: Residual analysis and diagnostics for GEE Models in R

text: Some colleagues asked me to perform a residual analysis on both linear models and generalized estimating equation (GEE) models. I know it is a faux-pas in some circles to remove outliers, but in our field it is common as a sensitivity analysis. In R, there are some great tools as diagnostics for linear models in the influence package using the `rstudent()` function. I am wondering if there are any suggestions on how to perform a similar studentized residual test for GEE models? I am also wondering the validity of doing this, considering the residuals are not independent? A previous poster asked a similar question and was directed to the `geeglm()` function in the geepack package. (That post is here: Diagnostics for GEE in R.) When I use `geeglm()` on my GEE model and try to use any of the functions in the influence package like `rstudent()`, I get the following error:               > test.gee <- geeglm(fsiq ~peabody+depress, id=id, data=df_long,                           family="gaussian", corstr="exchangeable")     > rstudent(test.gee, type="pearson")     Error in match.arg(type) :      'arg' should be one of “pearson”, “working”, “response”      Not sure if it's a coding error or error in my logic. Any help would be greatly appreciated.


title: Interpreting standardized mean centers in a cluster

text: I created a $k$-means with 3 clusters. Some of the variables had a big scale, so I used a $z$-score to standardize them. The others (mostly dummies), I left as is. Now, when I create the table of all variable centers for each cluster, obviously the ones I standardized are still expressed in the $z$-score format. The others are in their non-standardized state. So, how do I interpret the $z$-scored centers into a practical meaning (lets say I had normalized income, and its standardized center was -1.2). Do I convert that -1.2 back into something? Is it directly interpretable?


title: Coefficient / model averaging to control for exogenous circumstances in prediction

text: I'm interested in exploring statistical models (or modifications thereof) designed to handle a specific type of problem. Due to my ignorance of statistical terminology, I can only describe this type of problem by (contrived) examples: Suppose we're interested in estimating the likelihood that a given cell phone customer will drop his/her service (churn) over the course of the next month (October). We are given cell phone user data from the last three months (July, August, September). In estimating churn probability for October, we'll want to make use of a model that emphasizes inputs that are consistently predictive of churn each month, suppressing the importance of those whose coefficients vary heavily on a month-by-month basis. For example, "market" may be consistently predictive of churn, but an externality specific to one of the training datasets (e.g., marketing of new phone by a competitor in a specific market in August) may affect a coefficient in a way that won't necessarily apply to October. The market coefficient should be corrected/compressed/averaged when predicting October's churn likelihood. Here's a less contrived example. In estimating the likelihood that a voter will vote in the 2012 election, we may use voter-level data from the 2000, 2004, and 2008 elections to train a logistic regression on turnout. If "party" is an input, 2000 and 2004 coefficients may differ significantly from 2008 party coefficient due to exogenous circumstances (a unique political environment). In estimating 2012 vote likelihood for a given voter, we'll want to trim/compress/average the party coefficient to "smooth out" election-year- specific differences that can't be easily quantified. Can anyone point me in the right direction? Any guidance on how to better ask this question? What terminology should I be using? Thanks in advance for your help.


title: How can I perform a post hoc test of a three-way interaction in SPSS using the 'test' command?

text: I will describe the design of the experiment briefly.   1. The task have two within-subject variables ($A$ and $B$),  A. have two levels B. have three levels   2. one repeated variable with 5 levels ('time lag' 1,2,3,4,5).   3. Subjects are two groups. The problem is $A$ and $B$ were repeated with 'time lag' because its the condition of stimuli. i.e. A1B1lag3 A1B3lag5 .... Below is the SPSS syntax code.                   MIXED Accuracy BY Lag A B Group            /CRITERIA=CIN(95) MXITER(100) MXSTEP(100) SCORING(1) SINGULAR(0.000000000001) HCONVERGE(0,              ABSOLUTE) LCONVERGE(0, ABSOLUTE) PCONVERGE(0.000001, ABSOLUTE)            /FIXED=Lag A B Group Lag*A Lag*B Lag*Group A*B A*Group B*Group Lag*A*B Lag*A*Group Lag*B*Group A*B*Group Lag*A*B*Group           | SSTYPE(3)            /RANDOM=INTERCEPT  | SUBJECT(SubID) COVTYPE(ID)            /METHOD=REML          /repeated = **Lag*A*B** | subject(SubId) covtype(UN)          /PRINT = SOLUTION TESTCOV LMATRIX      The design works. However, when I try to use the /Test command to compare Estimated means, it always gave me the error "Custom hypothesis test 1 will not be performed because the L matrix is not estimable." The command /Test is something like below.                /test = 'Emotion 12v3' Emotion -.5 -.5 1         Agency*Emotion -1/4 -1/4 1/2 -1/4 -1/4 1/2         Agency*Emotion*Group -1/4 -1/4 1/2 -1/4 -1/4 1/2 1/4 1/4 -1/2 1/4 1/4 -1/2      I tried to test whether there is a group difference in A*B interaction… Also know that I can use EMMEANS               /EMMEANS=TABLES(Group*A*B) COMPARE ADJ(BONFERRONI)      However, I prefer 'test'... If I am wrong (or test and EMMEANS are eventually same one...) please advice me. Thanks.


title: Find out if there is a way to create BMI weight category variable

text: I am working on a project and have to find a way to break my data into percentiles. **The variables that I already have are gender, age in months and BMI. If I am able to calculate the percentile that each child falls into, I will then be able to label them as underweight, normal weight, overweight and obese.** I would be using the following cutoff values: > CDC Growth Charts 5th and 95th percentile > > BMI-for-age ≥ 95th Obesity   >  BMI-for-age ≥ 85th and < 95th Overweight   >  BMI-for-age < 5th Underweight The previous links weren't working, so here are links to the growth charts for boys (pdf) and girls (pdf) if anyone isn't familiar with them. These charts present contours of BMI percentiles (from the 5th through the 95th in irregular intervals) by age (for ages 24 to 240 months). Can anyone suggest a simple and accurate way to use these charts to obtain the preceding classifications from my data?


What is the meaning of $\|a\|_p=\left(\sum _{i=1}^n \left|a_i(t)\right|{}^p\right){}^{\frac{1}{p}}$? This formula is called out on the fifth page of _An Improved Data Stream Summary: The Count-Min Sketch and its Applications_ (which can be found here). I'm implementing the Count-Min Sketch and can understand the basic concepts just fine, but some of the finer points are explained in terms of this equation and some other terminology that I'm unacquainted with.


title: Logistic / multinomial regression as two / multiple Poisson regressions?

text: Can we instead of doing logistic or multinomial regression do two or multiple Poisson regressions and then combine Poisson predictions to get probabilistic predictions? If yes, how should we transform data for training?


title: How to properly analyse longitudinal data from experiments?

text: I am wondering how one should _properly_ analyse longitudinal data from simple experiments (not observational studies). Imagine an experiment where subjects are first measured at baseline, randomised to different treatments and then followed for a number of time periods. We are interested in differences (typically in means, but perhaps also variances) between the the effect of treatments at each time period. We only look at a single outcome measure. # The simple(st) case Let $Y_{ij}$ be the measurement of subject _i_ at time period _j_ , $i=1,\dotsc,n$, $j=0,\dotsc,k$, where $j=0$ corresponds to baseline (i.e., pre-intervention) measurements. If we have two treatments and only $k=1$ follow-up time period (and are interested in means), it’s easy to analyse the data. We can either look at the change scores $Y_{i1}-Y_{i0}$ and compare these between the two treatment groups using a _t_ -test, or do an ANCOVA where we adjust for baseline, $Y_{i1}=\alpha + \beta_1 Y_{i0} +\beta_2 T_i+\epsilon_i$, where $T_i$ is either 0 or 1, depending on the treatment subject _i_ received. Both methods of analysis are valid, but the ANCOVA is more powerful (and may give us ‘better’ answers if the treatment groups happened to be unbalanced at baseline?). (Note that all this for an _experiment_ ; for observational studies, other analyses may be more appropriate.) # The more interesting/difficult case But how should we properly analyse data like this when we have measurements at _more_ than one follow-up time? It seems to be that things quickly gets much more complicated. The simplest example is an experiment with _two_ follow-up times. Let’s make it more concrete, and imagine that we are looking at weight- loss, and are investigating the effect of an intervention (diet, exercise or nutritional advice) on overweight people. So we gather a bunch of overweight people, measure their weight (or BMI), randomise them to either intervention or control, perform the intervention, and then measure them again after, say, 6 months and 2 years. I see several ways to analyse this:   1. Ignore the baseline values and analyse only the follow up-data.   2. Compare change scores like we did in the simpler case. We now have two changes scores $Y_{i1}-Y_{i0}$ and $Y_{i2}-Y_{i0}$.   3. Model the three measurements as correlated random variables, and estimate the parameters using generalised least squares with an unstructured covariance matrix.   4. Do an ANCOVA-style analysis at each time point. Here’s my understanding of these methods:   1. This seems like a really bad idea, as we’re throwing away valuable data. :-( But it _is_ a valid method of analysis.   2. This also seems valid, but I guess it too wastes power. And we really need to jointly model to two change scores. How should we do this?   3. This method of analysis initially sounded perfect (to me). We get nice estimates for both the means at each time point and the difference in means between the two groups at each time point. We can easily add the restriction that the mean at baseline is equal for the two groups (to gain some power), but one problem is the distributional assumptions. Even if we are prepared to assume approximately normally distributed errors at the two follow-up times, for the baseline measurements the errors are guaranteed _not_ to be normally distributed, as the inclusion criterion was that the people participating should be overweight (e.g., weight > some fixed constant). So if weight was normally distributed in the population, it would have a _truncated_ normal distribution for the baseline measurements for our subjects (at follow-up, this would probably be less of a problem). A similar issue would of course arise if the inclusion criteria didn’t _exactly_ correspond to limits on baseline values, but only were (cor)related with them (e.g., where the inclusion criteria specifies a BMI > a constant, and we are measuring triglycerides). Because of this, it doesn’t feel right treating the baseline measurements on the same footing as the follow-up measurements; they are fundamentally different.   4. This I’m actually not sure how one should go about doing. I guess we could create a linear model $Y_{ij}=\beta_j + \beta_{1j} Y_{i0} +\beta_{2j} T_{i} + \epsilon_{ij}$, $j=1, 2$, and fit this using generalised least squares (with an unstructured covariance matrix). Are there any conceptual or technical problems with doing this? Will a similar method work just as well even if we have many different treatments and/or time points (as long as we have plenty of observations)? Or if we have some missing data (MCAR or MAR)? Basically, my question is: For longitudinal data in experiments like the one described above, which of the methods 1, 2, 3 and 4 (or other methods) can be considered _valid_ methods of analyses, and which of them are _preferred_ (e.g., have the greatest power, minimum variance of estimators, least bias …)? I know many different methods of analysis have been used in the past, but what are the current recommendations?


title: Where to cut a dendrogram?

text: Hierarchical clustering can be represented by a dendrogram. Cutting a dendrogram at a certain level gives a set of clusters. Cutting at another level gives another set of clusters. How would you pick where to cut the dendrogram? Is there something we could consider an optimal point? If I look at a dendrogram across time as it changes, should I cut at the same point?


title: Feature selection methods for document classtification

text: I have a simple document classification problem where i need to classify some documents to a definite set of classes. I need to perform a feature selection (where I will select the most important words from each class). Currently I am calculating the tf-idf scores. Is there a better way to do that ? I heard about chi-square statistic being used in this context, is that true? And if that is the case, could you give me (preferably simple) links to resources that I can read for further information about this? thanks


title: When/where to use functional data analysis?

text: I am _very_ new to functional data analysis (FDA). I am reading: > Ramsay, James O., and Silverman, Bernard W. (2006), Functional Data > Analysis, 2nd ed., Springer, New York. However, I am still not very clear where/when to use FDA? Could someone please give me an example especially in medical studies? I really don't know where/when to apply FDA in practice. For growth curve data, we can use nonlinear mixed models, for longitudinal data, we can use repeated measure ANOVA, and for multivariate data/high dimensional data, we can use PCA, FA, etc. So when/where will be the best timing/situation to use FDA?


title: Statistics for small and unequal samples

text: I have a pre-test, post-test, non-equivalent groups experiment, with a sample of 62 (control and experimental groups have 31 each). My objective is to compare the effectiveness of a model of teaching for developing vocabulary among pupils with different learning styles. The learning styles of the 31 of the experimental group are:   * Visual = 5,    * Auditory = 7,    * Reading = 12,    * Kinaesthetic = 1,    * VA = 3,    * VR = 1,    * AR = 2.  What statistics should I employ to compare their performance (scores on vocabulary test)?


title: How do I plot an exponential curve in R?

text: I am trying to plot an exponential curve (`nls`) through this data set in R. abm is a text file with the following data=               x=          1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,2,1,2,2,1,2,1,2,2,2,2,2,1,2,2,1,3,1,3,1,2,1,2,2,3,3,2,3,3,3,2,3,3,2,1,5,5,1,4,3,4,3,3,5,7,2,4,6,4,4,3,6,5,5,3,4,3,4,7,5,7,11,3,6,6,6,9,1,1,3,9,3,4,6,6,12,6,2,12,8,2,4,12,6,8,9,6,5,4,4,11,4,18,8,10,8,17,5,16,16,6,5,18,3,1,10,2,11,21,8,11,15,9,13,14,3,11,17,18,9,16,17,12,18,2,9,17,14,14,5,21,21,13,17,15,9,1,18,17,7,2,17,23,22,27,10,10,17,25,12,28,12,13,14,5,20,13,3,17,16,21,11,27,19,19,19,22,19,6,30,16,13,10,30,18,17,30,4,23,20,18,5,22,10,37,38,19,36,36,34,37,36,25,30,26,22,32,30,30,45,18,25,42,44,24,38,34,14,16,6,36,26,29,28,44,20,31,45,39,43,46,47,34,46,47,40,48,43,49,46,48,49,49          y=          18.36734694,18.36734694,18.36734694,18.36734694,18.36734694,18.36734694,18.36734694,18.36734694,18.36734694,18.36734694,20.40816327,20.40816327,20.40816327,20.40816327,20.40816327,20.40816327,22.44897959,22.44897959,22.44897959,22.44897959,24.48979592,24.48979592,26.53061224,26.53061224,26.53061224,28.57142857,28.57142857,30.6122449,30.6122449,32.65306122,32.65306122,34.69387755,34.69387755,36.73469388,38.7755102,38.7755102,38.7755102,40.81632653,40.81632653,40.81632653,40.81632653,40.81632653,42.85714286,44.89795918,48.97959184,51.02040816,51.02040816,53.06122449,55.10204082,57.14285714,61.2244898,61.2244898,63.26530612,63.26530612,63.26530612,67.34693878,67.34693878,69.3877551,77.55102041,77.55102041,77.55102041,79.59183673,81.63265306,81.63265306,89.79591837,89.79591837,95.91836735,104.0816327,112.244898,124.4897959,130.6122449,134.6938776,138.7755102,146.9387755,146.9387755,151.0204082,161.2244898,163.2653061,165.3061224,167.3469388,169.3877551,169.3877551,169.3877551,171.4285714,181.6326531,185.7142857,208.1632653,220.4081633,224.4897959,230.6122449,240.8163265,240.8163265,246.9387755,267.3469388,269.3877551,271.4285714,304.0816327,306.122449,310.2040816,318.3673469,332.6530612,342.8571429,355.1020408,355.1020408,387.755102,395.9183673,395.9183673,397.9591837,404.0816327,430.6122449,448.9795918,451.0204082,455.1020408,465.3061224,465.3061224,467.3469388,469.3877551,518.3673469,522.4489796,528.5714286,532.6530612,534.6938776,534.6938776,553.0612245,563.2653061,569.3877551,577.5510204,583.6734694,583.6734694,600,612.244898,618.3673469,628.5714286,689.7959184,714.2857143,730.6122449,738.7755102,748.9795918,802.0408163,855.1020408,875.5102041,895.9183673,906.122449,924.4897959,934.6938776,948.9795918,955.1020408,977.5510204,989.7959184,997.9591837,1014.285714,1026.530612,1040.816327,1071.428571,1083.673469,1102.040816,1110.204082,1126.530612,1130.612245,1142.857143,1185.714286,1218.367347,1224.489796,1236.734694,1242.857143,1248.979592,1259.183673,1314.285714,1424.489796,1430.612245,1467.346939,1495.918367,1502.040816,1528.571429,1534.693878,1540.816327,1557.142857,1608.163265,1632.653061,1657.142857,1691.836735,1691.836735,1738.77551,1765.306122,1916.326531,1928.571429,2055.102041,2057.142857,2059.183673,2071.428571,2104.081633,2155.102041,2171.428571,2193.877551,2226.530612,2251.020408,2357.142857,2508.163265,2553.061224,2579.591837,2581.632653,2738.77551,2740.816327,2759.183673,2897.959184,2902.040816,2904.081633,2914.285714,3057.142857,3193.877551,3204.081633,3306.122449,3373.469388,3381.632653,3428.571429,3491.836735,3604.081633,3875.510204,4028.571429,4273.469388,4781.632653,5002.040816,5004.081633,5575.510204,5683.673469,6073.469388,6185.714286,6212.244898,6391.836735,6434.693878,6442.857143,6591.836735,6591.836735,6724.489796,6775.510204,8030.612245,8710.204082,9410.204082,9434.693878,9991.836735,10787.7551,11000,11048.97959,11524.4898,12085.71429,12881.63265,13620.40816,14087.7551,14089.79592,15636.73469,16361.22449,16440.81633,16759.18367,16846.93878,18267.34694,20240.81633,22004.08163,22816.32653,26493.87755,29981.63265,32777.55102,33724.4898,37240.81633,39116.32653,54540.81633,56977.55102,73912.2449      I use this R script to make a scatter plot:               plot(log(abm), xlab="Log10 (Number of sites occupied)", ylab="(Log10) Mean local abundance", xlim=c(0,4),pch=20)      Which looks like this: ![enter image description here](http://i.stack.imgur.com/c8RcM.jpg) Now I want to plot an exponential curve through this data. Can anybody please help with this? I know I need to use "nls" but I cannot seem to do it. I am a noob at R and would appreciate any advice and help. Cheers!


title: "Array index is less than one" WINBUGS Error

text: for (i in 1:160) {         situps[i] ~ dnorm(mu[i],prec[gender[i]+1]);         mu[i] <- b0+b1[gender[i]]+b2[health[i]]+b1b2[health[i],gender[i]];     }     b0 ~ dnorm(0,.001);     for (i in 1:4) {         b2[i] ~ dnorm(0,.001);     }     for (j in 1:2) {         b1[j] ~ dnorm(0,.001);         prec[j] <- pow(sigma[j],-2);         sigma[j] ~ dunif(0,10);     }     for (i in 1:4) {         for (j in 1:2) {             b1b2[i,j] ~ dnorm(0,.001);         }     }      I'm not sure why this gives me the error "array index is less than one" because I clearly have all of my arrays starting at 1. Any ideas?


title: Programmer looking to break into machine learning field

text: I am a software developer (mostly .NET and Python about 5 years experience). What can I do to help me get a job in the machine learning field or really anything that will get me started in that field? Is post-graduate degree a hard requirement?