How to find the variability of your estimate (aka standard error)?

summary of the lm function will give u the answer, but it's standard error of the TRAINING data.

for X1, that's 0.02593

How to find the standard error of TEST data?

bootstrap - concept is to

1. randomly sample the TRAINING data with replacement to produce a new data set
2. do the above over and over again and you'll have many data sets
3. for each data set, you get an estimate and you'll end up with many estimates
4. standard deviation of the sample estimates will be your standard errors

(quite clever!)

> summary(lm(y~X1+X2, Xy))

Call:
lm(formula = y ~ X1 + X2, data = Xy)

Residuals:
     Min       1Q   Median       3Q      Max 
-1.44171 -0.25468 -0.01736  0.33081  1.45860 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  0.26583    0.01988  13.372  < 2e-16 ***
X1           0.14533    0.02593   5.604 2.71e-08 ***
X2           0.31337    0.02923  10.722  < 2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 

Residual standard error: 0.5451 on 997 degrees of freedom
Multiple R-squared: 0.1171, Adjusted R-squared: 0.1154 
F-statistic: 66.14 on 2 and 997 DF,  p-value: < 2.2e-16 

> coef(summary(lm(y~X1+X2, Xy[1:10,])))
             Estimate Std. Error   t value     Pr(>|t|)
(Intercept) -5.691507  1.4521776 -3.919292 0.0057542755
X1          -3.573758  0.6417983 -5.568350 0.0008434796
X2          13.189553  2.8350903  4.652252 0.0023356811
> coef(summary(lm(y~X1+X2, Xy[1:10,])))["X1", "Estimate"]
[1] -3.573758

> est.fn=function(data, index){
+ coef(summary(lm(y~X1+X2, data[index,])))["X1", "Estimate"]}
> est.fn(Xy, 1:10)
[1] -3.573758
> est.fn(Xy, 1:1000)
[1] 0.1453263
> boot.out=boot(Xy,est.fn,R=1000)
> boot.out

ORDINARY NONPARAMETRIC BOOTSTRAP


Call:
boot(data = Xy, statistic = est.fn, R = 1000)


Bootstrap Statistics :
     original       bias    std. error
t1* 0.1453263 0.0002518276   0.0302769
> plot(boot.out) 

for time series, ie the data are not idd, we do boostrap in a block of certain size:

> tsboot(Xy, est.fn, 1000, sim="fixed", l=100)
BLOCK BOOTSTRAP FOR TIME SERIES
Fixed Block Length of 100 
Call:tsboot(tseries = Xy, statistic = est.fn, R = 1000, l = 100, sim = "fixed")

Bootstrap Statistics :     original      bias    std. errort1* 0.1453263 0.001404265   0.1983411

How precise are your regression estimators?

Short answer: standard error

In details:
Let's say we wanna estimate the population mean from sample data.
It makes sense to expect the sample mean to be very close to population mean.
Furthermore, if we keep getting more new sample sets, we would expect the error of our estimate to be smaller.
Taking it one step further, if the means from those sample data set don't differ much, we'd expect the error of estimate to be smaller and vice versa.

So, roughly speaking, standard error of the mean estimate is
* proportional to how much each sample  value deviates from sample mean (actually, square root of sum of square difference)
* inversely proportional to sample size (actually,  square root of n-2)

How about the standard error of slope estimate in linear regression?

roughly speaking, standard error of the slope estimate is
* proportional to how much each regression predicted value deviates from observed value  (actually, square root of sum of square difference)
* inversely proportional to sum of square difference of PREDICTOR.
Another words, wider range of your predictor, smaller your standard error .
It actually makes sense. Think about it. Slope = (y2-y1)/(x2-x1)
If y2 and/or y1 is off a little bit
* how much the slope is gonna change given x2 is far away from x1?
* how much the slope is gonna change given x2 is very close to x1?
So, if u are doing a controlled experiment, choose a wider range of predictors.

Closely related is t-statistics,  which is defined as
Estimated value/ std error
With a given t-statistic, we can look up a corresponding p-value, which is the probability of obtaining the data if the null hypothesis is true.
95% confidence interval is about estimated value +- 2*std error
If the interval contains 0, you cannot reject the null hypothesis.

How to determine how good your regression model is?

It comes down to model variance - model bias trade off.

What's model variance?
Let's have our model fit the data as much as possible.
If we randomly pick a certain % of our data and fit as much as possible,  we have model 1.
We do this over and over again for n times and end up with n models.
The model variance is basically the variance of the predicted value (y) of these n models for a given x.
Typically,  if we use a less flexible model, say linear, the model variance is gonna be less.

What's model bias?
That's the difference between average predicted y from n models for a given x and the actual y for the same x.
Typically,  if we use a less flexible model, say linear, the model bias is gonna be bigger.

We'd like to arrive at a better model by trading off between model variance and model bias. We don't want our model be too flexible to overfit the data and introduce a huge model variance while we don't want our model be to rigid and have the predicted value be too far away from actual value.

How to tell if u are overfitting?
Use out of sample data.
Let's say u fit the training data and come up with a model with a very small mean sq error on the training data.
Have this model to predict the out of sample data and compare with the actual result.  If u have a large mean sq error on the out of sample data, you are probably overfitting your training data.

MySQL Stored Procedures

visit: http://dev.mysql.com/tech-resources/articles/mysql-storedprocedures.pdf

Linux VM + MySQL + R + DBI

• Linux VM 
• get it @c9.io

• MySQL
• type mysql on c9.io terminal and you'll have mysql server access

• R
• if old version is sufficient (note, DBI requires >= 2.51), run
• c9pm install r-lang
• or, manual install for newer version (take quite a while...), see:
• http://stackoverflow.com/questions/13146852/r-2-15-install-in-redhat
• http://cran.r-project.org/doc/manuals/r-release/R-admin.html#Installing-R-under-Unix_002dalikes

• DBI
• from R, install.packages('DBI')
• RMySQL
• wget 'http://cran.r-project.org/src/contrib/RMySQL_0.9-3.tar.gz'
• R CMD INSTALL RMySQL_0.9-3.tar.gz 

Try it from R:

> library(DBI)

> library(RMySQL)

> con <- dbConnect(MySQL(), user="USER", dbname="DB", host="HOSTNAME")

> dbGetQuery(con, "select curdate() from dual;")

   curdate()

1 2014-01-03

Null hypothesis test and Type I / Type II errors

Null hypothesis test

- why do it?  to have some idea whether there is a relationship between an independent variable (IV) and a dependent variable (DV)

- how? 

1. make the hypothesis that there is NO relationship (null hypothesis)
2. calculate the p-value, which is the PROBABILITY OF GETTING THESE DATA GIVEN NULL HYPOTHESIS IS TRUE, rather than the p(null hypothesis is true)
3. If p-value is small enough (by convention under 5%, which is called α (alpha) level), we'll reject the null and claim that there's a statistical significant relationship between the IV and DV

Type I error

- false positive
- how to remember?

1. positive <- hey, there's a relationship!
2. false <- actually there is not... :(
3. that means we have rejected the null while we shouldn't have 

- α (alpha) = P(Type I error | null hypothesis is true)
 

Type II error

- false negative
- i.e. failed to reject the null 
- β (beta) = P(Type II error | null hypothesis is false)

Power

= P(correctly reject the null) = 1- β


Ref: http://onlinestatbook.com/2/logic_of_hypothesis_testing/errors.html