df=read.table("d:/lin/teaching/NDHU/Econometric Analysis/data.txt",header=TRUE)
df$Y
df[,2]
lmout=lm(Y ~ X1 + X2, data=df)
names(lmout)
print(lmout)
summary(lmout)
lmout$co
str(df)
 df$X1
 length(df$X1)
 is.vector(df$X1)
 vec=c(1,2,3,4,5,6)
 vec
 is.vector(vec)
 mat=matrix(c(1,2,3,4,5,6),ncol=2)
 mat
 mat=matrix(c(1,2,3,4,5,6),byrow=TRUE,ncol=2)
 mat
 dfmat=as.matrix(df)
 dfmat
 dim(dfmat)
 mat1=matrix(c(1,2,3,4,5,6),byrow=TRUE,ncol=2)
 mat2=matrix(c(7,8,9,10,11,12),,ncol=2)
 mat1
 mat2
 mat3=cbind(mat1,mat2)
 mat4=rbind(mat1,mat2)
 mat3
 mat4
 rbind(mat1,c(99,99))
 ar=array(c(1,2,3,4,5,6),dim=c(3,2,2))
 ar
 l=list(1,"a",c(4,4),list(FALSE,2))
 Lst <- list(name="Fred", wife="Mary", no.children=3,child.ages=c(4,7,9))
 
 # Accessing Elements and Subsetting Vectors, Arrays, and Lists
 vec=c(1,2,3,2,5)
 vec[3]
 index=c(3:5)
 index
 vec[index]
 index=vec==2
 index
 vec[index]

vec[vec!=2]
vec[-c(3:5)]

dfmat
ldata=list(A=dfmat[dfmat[,1]=="A",2:4],B=dfmat[dfmat[,1]=="B",2:4])
ldata
ldatadf=list(A=df[df[,1]=="A",2:4],B=df[df[,1]=="B",2:4])

 ldatadf=list(A=df[df[,1]=="A",2:4],B=df[df[,1]=="B",2:4])
 lmout=NULL
 for (i in 1:2) {
 lmout[[i]]=lm(Y ~ X1+X2,data=ldatadf[[i]])
 print(lmout[[i]])
 }
 
apply(df[,2:4],2,mean)
 
y=as.numeric(dfmat[,2])
X=matrix(as.numeric(dfmat[,3:4]),ncol=2)
X=cbind(rep(1,nrow(X)),X)
XpXinv=chol2inv(chol(crossprod(X)))
bhat=XpXinv%*%crossprod(X,y)
res=as.vector(y-X%*%bhat)
ssq=as.numeric(res%*%res/(nrow(X)-ncol(X)))
se=sqrt(diag(ssq*XpXinv))
 
myreg=function(y,X){
#
# purpose: compute lsq regression
#
# arguments:
# y -- vector of dep var
# X -- array of indep vars
#
# output:
# list containing lsq coef and std errors
#
XpXinv=chol2inv(chol(crossprod(X)))
bhat=XpXinv%*%crossprod(X,y)
res=as.vector(y-X%*%bhat)
ssq=as.numeric(res%*%res/(nrow(X)-ncol(X)))
se=sqrt(diag(ssq*XpXinv))
list(b=bhat,std_errors=se)
}
 myreg(X=X,y=y)
 
 #
# purpose: compute lsq regression
#
# arguments:
# y -- vector of dep var
# X -- array of indep vars
#
# output:
# list containing lsq coef and std errors
#
XpXinv=chol2inv(chol(crossprod(X)))
bhat=XpXinv%*%crossprod(X,y)
res=as.vector(y-X%*%bhat)
ssq=as.numeric(res%*%res/(nrow(X)-ncol(X)))
se=sqrt(diag(ssq*XpXinv))
cat("in myreg, se = ",fill=TRUE)
print(se)
list(b=bhat,std_errors=se)
}

 myregout=myreg(y,X)

 debug(myreg)
 myreg(X,y)
 undebug(myreg)
 
 # Graphics
 hist(rnorm(1000),breaks=50,col=ˇ¨magentaˇ¨)
 hist(rnorm(1000),breaks=30,col=ˇ¨magnetaˇ¨,xlab=ˇ¨thetaˇ¨,ylab=ˇ¨ˇ¨,
 main=ˇ¨Non-parametric Estimate of Theta Distributionˇ¨)

 par(mfrow=c(2,2))
X=matrix(rnorm(5000),ncol=5)
X=t(t(X)+c(1,4,6,8,10))
hist(X[,1],main="Histogram of 1st col",col="magenta",xlab="")
plot(X[,1],X[,2],xlab="col 1", ylab="col 2",pch=17,col="red",
xlim=c(-4,4),ylim=c(0,8))
title("Scatterplot")
abline(c(0,1),lwd=2,lty=2)
matplot(X,type="l",ylab="",main="MATPLOT")
acf(X[,5],ylab="",main="ACF of 5th Col")

memory.size()
 getwd()
 x=matrix(rnorm(1e07),ncol=1000)
 memory.size()
 memory.limit()
 begin=proc.time()[3]
 z=crossprod(x)
 end=proc.time()[3]
 print(end-begin)

test=function(n){x=matrix(rnorm(n),ncol=1000);z=crossprod(x);
cz=chol(z)}
 Rprof("test.out")
 test(1e07)
 Rprof()
 summaryRprof("test.out")
 
 n=1e04
 x=NULL
 zero=c(rep(0,5))
 begin=proc.time()[3]
 for (i in 1:n) {x=rbind(x,zero) }
 end=proc.time()[3]
 print(end-begin)
 x=NULL
 begin=proc.time()[3]
 x=matrix(double(5*n),ncol=5)
 end=proc.time()[3]
 print(end-begin)