rm(list=ls())
library(MASS)
m=c(1,1,1,1,1)
sigma=diag(m)
alpha=20
gamma=c(1,-0.5,0.7,2.5,-0.8)



f=function(n,J,K,l,t){

zz=matrix(NA,n,5);
for(i in 1:n){
zz[i,]=mvrnorm(1, rep(0, 5), sigma)}

u1=matrix(NA,J,n);
for(i in 1:n){
 for(j in 1:J){
u1[,i]=rnorm(J,0,(10/(j-0.5)*pi)^2)}}

u2=matrix(NA,J,n);
for(r in 1:n){
  for(j in 1:J){
u2[,r]=rnorm(J,0,(10/(j-0.5)*pi)^2)}}


u=c();
ss=matrix(NA,n,K);
xt=array(NA,dim=c(J,K,n))
h=runif(K,0,1)
x=c();x1=c();
X=function(K,J){
 for(i in 1:n){
  for(k in 1:K){
   for(j in 1:J){
     xt[j,k,i]=sqrt(2)*(u1[,i])[j]*sin((j-0.5)*pi*h[k])}}}

 for(l in 1:n){
u=c(u,apply(xt[,,l],2,sum))}
ss=matrix(u,n,K)
ss}
XX=X(K,J)

b=c();
beta=function(K){
 for(k in 1:K){
  b[k]=sqrt(2)*sin(pi*h[k]/2)+2*sqrt(2)*sin(3*pi*h[k]/2)}
b}
bb=beta(K)


xb=matrix(NA,n,K);
XB=function(K,J,n){
 for(i in 1:n){
  for(k in 1:K){
    xb[i,k]=(XX[i,k])*(bb[k])}}
xb}

riman=c();
for(i in 1:n){
riman[i]=sum(XB(K,J,n)[i,])/K}


e=rnorm(n)

y=c();
Y=function(n){
for(i in 1:n){
 y[i]=alpha+t(zz[i,])%*%as.matrix(gamma)+riman[i]+e[i]}
y}


x3=c();
Sigma=function(n,J){
s=matrix(NA,n,n)
for(r in 1:n){
for(i in 1:n){
 for(j in 1:J){
x3[j]=((u1[,i])[j])*((u2[,r])[j])*((-2*cos(j-0.5)*pi/(j-0.5)*pi)+(2*sin(j-0.5)*pi/((j-0.5)*pi)^2)
-((4/3)*(1/((j-0.5)*pi)^2))*((cos(j-0.5)*pi)-1))}
s[i,r]=sum(x3)}}
s}

x1=c();t1=c();
T1=function(n,J){
for(i in 1:n){
for(j in 1:J){
x1[j]=(u1[,i])[j]*sqrt(2)*(1/(j-0.5)*pi)*(1-cos((j-0.5)*pi))}
t1[i]=sum(x1)}
t1}

x2=c();y1=c();t2=c();
T2=function(n,J){
for(i in 1:n){
for(j in 1:J){
     x2[j]=-(u1[,i])[j]*sqrt(2)*(1/(j-0.5)*pi)*cos((j-0.5)*pi)
     y1[j]=(u1[,i])[j]*sqrt(2)*(1/((j-0.5)*pi))^2*sin((j-0.5)*pi)}
t2[i]=(sum(x2)+sum(y1))}
t2}


T=matrix(NA,n,8)
for(i in 1:n){
T[i,1]=1
T[i,2]=as.matrix(T1(n,J))[i,1]
T[i,3]=as.matrix(T2(n,J))[i,1]
T[i,4:8]=zz[i,]
}


lambda=c(10^-10,10^-9,10^-8,10^-7,10^-6,10^-5,10^-4,10^-3,10^-2,10^-1,1,10,10^2,10^3)
M=function(n,l){
   a=rep(1,n)
   diag(a)
   m=Sigma(n,J)+n*lambda[l]*diag(a)
m}


d1=solve(t(T)%*%solve(M(n,l))%*%T)%*%t(T)%*%solve(M(n,l))%*%as.matrix(Y(n))

c1=solve(M(n,l))%*%(1-T%*%solve(t(T)%*%solve(M(n,l))%*%T)%*%t(T)%*%solve(M(n,l)))%*%as.matrix(Y(n))

d=t(d1)
c=t(c1)


x4=c();y2=c();x5=c();xk=c();
XK=function(t){
for(i in 1:n){
 for(j in 1:J){
   x4[j]=-(u1[,i])[j]*sqrt(2)*(1/(j-0.5)*pi)*cos((j-0.5)*pi)
   y2[j]=(u1[,i])[j]*sqrt(2)*(1/((j-0.5)*pi))^2*sin((j-0.5)*pi)
   x5[j]=(u1[,i])[j]*sqrt(2)*(1/(j-0.5)*pi)*(1-cos((j-0.5)*pi))}
xk[i]=(t-1/2)*(sum(x4)+sum(y2))+(1/3-t/2)*sum(x5)}
xk}


betahat=function(t){
beta=d[2]+d[3]*t+sum(c*XK(t))
beta}


beta=function(t){
sqrt(2)*sin(pi*t/2)+2*sqrt(2)*sin(3*pi*t/2)}

##############yhat
bbetahat=function(K){
for(k in 1:K){
bhat=d[2]+d[3]*t+sum(c*XK(K))}
bhat}
bbhat=bbetahat(K)


xbhat=matrix(NA,n,K);
XBhat=function(K,J,n){
 for(i in 1:n){
  for(k in 1:K){
    xbhat[i,k]=(XX[i,k])*(bbhat[k])}}
xbhat}

rimanhat=c();
for(i in 1:n){
rimanhat[i]=sum(XBhat(K,J,n)[i,])/K}

yhat=c();
Yhat=function(n){
for(i in 1:n){
 yhat[i]=d[1]+t(zz[i,])%*%as.matrix(d[4:8])+rimanhat[i]}
yhat}



#################

c(betahat(t),d[1],d[4:8],beta(t))
}


#########################100

biasbeta=mean(replicate(1000,f(100,100,100000,14,0.2)[1]))-f(100,100,100000,14,0.2)[2]
varbeta=var(replicate(1000,f(100,100,100000,14,0.2)[1]))
msebeta1=varbeta+biasbeta^2



biasalpha=mean(replicate(1000,f(100,100,100000,14,0.2)[2]))-alpha
varalpha=var(replicate(1000,f(100,100,100000,14,0.2)[2]))
msealpha1=varalpha+biasalpha^2



biasgamma1=mean(replicate(1000,f(100,100,100000,14,0.2)[4]))-gamma[1]
vargamma1=var(replicate(1000,f(100,100,100000,14,0.2)[4]))
msegamma11=vargamma1+(biasgamma1)^2

biasgamma2=mean(replicate(1000,f(100,100,100000,14,0.2)[5]))-gamma[2]
vargamma2=var(replicate(1000,f(100,100,100000,14,0.2)[5]))
msegamma21=vargamma2+(biasgamma2)^2

biasgamma3=mean(replicate(1000,f(100,100,100000,14,0.2)[6]))-gamma[3]
vargamma3=var(replicate(1000,f(100,100,100000,14,0.2)[6]))
msegamma31=vargamma3+(biasgamma3)^2

biasgamma4=mean(replicate(1000,f(100,100,100000,14,0.2)[7]))-gamma[4]
vargamma4=var(replicate(1000,f(100,100,100000,14,0.2)[7]))
msegamma41=vargamma4+(biasgamma4)^2

biasgamma5=mean(replicate(1000,f(100,100,100000,14,0.2)[8]))-gamma[5]
vargamma5=var(replicate(1000,f(100,100,100000,14,0.2)[8]))
msegamma51=vargamma5+(biasgamma5)^2



################################200

biasbeta=mean(replicate(1000,f(200,100,100000,14,0.2)[1]))-f(200,100,100000,14,0.2)[2]
varbeta=var(replicate(1000,f(200,100,100000,14,0.2)[1]))
msebeta2=varbeta+biasbeta^2



biasalpha=mean(replicate(1000,f(200,100,100000,14,0.2)[2]))-alpha
varalpha=var(replicate(1000,f(200,100,100000,14,0.2)[2]))
msealpha2=varalpha+biasalpha^2



biasgamma1=mean(replicate(1000,f(200,100,100000,14,0.2)[4]))-gamma[1]
vargamma1=var(replicate(1000,f(200,100,100000,14,0.2)[4]))
msegamma12=vargamma1+(biasgamma1)^2

biasgamma2=mean(replicate(1000,f(200,100,100000,14,0.2)[5]))-gamma[2]
vargamma2=var(replicate(1000,f(200,100,100000,14,0.2)[5]))
msegamma22=vargamma2+(biasgamma2)^2

biasgamma3=mean(replicate(1000,f(200,100,100000,14,0.2)[6]))-gamma[3]
vargamma3=var(replicate(1000,f(200,100,100000,14,0.2)[6]))
msegamma32=vargamma3+(biasgamma3)^2

biasgamma4=mean(replicate(1000,f(200,100,100000,14,0.2)[7]))-gamma[4]
vargamma4=var(replicate(1000,f(200,100,100000,14,0.2)[7]))
msegamma42=vargamma4+(biasgamma4)^2

biasgamma5=mean(replicate(1000,f(200,100,100000,14,0.2)[8]))-gamma[5]
vargamma5=var(replicate(1000,f(200,100,100000,14,0.2)[8]))
msegamma52=vargamma5+(biasgamma5)^2




###################################400
biasbeta=mean(replicate(1000,f(400,100,100000,14,0.2)[1]))-f(400,100,100000,14,0.2)[2]
varbeta=var(replicate(1000,f(400,100,100000,14,0.2)[1]))
msebeta4=varbeta+biasbeta^2



biasalpha=mean(replicate(1000,f(400,100,100000,14,0.2)[2]))-alpha
varalpha=var(replicate(1000,f(400,100,100000,14,0.2)[2]))
msealpha4=varalpha+biasalpha^2



biasgamma1=mean(replicate(1000,f(400,100,100000,14,0.2)[4]))-gamma[1]
vargamma1=var(replicate(1000,f(400,100,100000,14,0.2)[4]))
msegamma14=vargamma1+(biasgamma1)^2

biasgamma2=mean(replicate(1000,f(400,100,100000,14,0.2)[5]))-gamma[2]
vargamma2=var(replicate(1000,f(400,100,100000,14,0.2)[5]))
msegamma24=vargamma2+(biasgamma2)^2

biasgamma3=mean(replicate(1000,f(400,100,100000,14,0.2)[6]))-gamma[3]
vargamma3=var(replicate(1000,f(400,100,100000,14,0.2)[6]))
msegamma34=vargamma3+(biasgamma3)^2

biasgamma4=mean(replicate(1000,f(400,100,100000,14,0.2)[7]))-gamma[4]
vargamma4=var(replicate(1000,f(400,100,100000,14,0.2)[7]))
msegamma44=vargamma4+(biasgamma4)^2

biasgamma5=mean(replicate(1000,f(400,100,100000,14,0.2)[8]))-gamma[5]
vargamma5=var(replicate(1000,f(400,100,100000,14,0.2)[8]))
msegamma54=vargamma5+(biasgamma5)^2

















<pre>An Error occurred while handling another error:
yii\web\HeadersAlreadySentException: Headers already sent in  on line 0. in /var/www/html/prof-homepages/vendor/yiisoft/yii2/web/Response.php:366
Stack trace:
#0 /var/www/html/prof-homepages/vendor/yiisoft/yii2/web/Response.php(339): yii\web\Response-&gt;sendHeaders()
#1 /var/www/html/prof-homepages/vendor/yiisoft/yii2/web/ErrorHandler.php(136): yii\web\Response-&gt;send()
#2 /var/www/html/prof-homepages/vendor/yiisoft/yii2/base/ErrorHandler.php(135): yii\web\ErrorHandler-&gt;renderException()
#3 [internal function]: yii\base\ErrorHandler-&gt;handleException()
#4 {main}
Previous exception:
yii\web\HeadersAlreadySentException: Headers already sent in  on line 0. in /var/www/html/prof-homepages/vendor/yiisoft/yii2/web/Response.php:366
Stack trace:
#0 /var/www/html/prof-homepages/vendor/yiisoft/yii2/web/Response.php(339): yii\web\Response-&gt;sendHeaders()
#1 /var/www/html/prof-homepages/vendor/yiisoft/yii2/base/Application.php(656): yii\web\Response-&gt;send()
#2 /var/www/html/prof-homepages/vendor/faravaghi/yii2-filemanager/models/Files.php(696): yii\base\Application-&gt;end()
#3 /var/www/html/prof-homepages/vendor/faravaghi/yii2-filemanager/controllers/FilesController.php(484): faravaghi\filemanager\models\Files-&gt;getFile()
#4 [internal function]: faravaghi\filemanager\controllers\FilesController-&gt;actionGetFile()
#5 /var/www/html/prof-homepages/vendor/yiisoft/yii2/base/InlineAction.php(57): call_user_func_array()
#6 /var/www/html/prof-homepages/vendor/yiisoft/yii2/base/Controller.php(180): yii\base\InlineAction-&gt;runWithParams()
#7 /var/www/html/prof-homepages/vendor/yiisoft/yii2/base/Module.php(528): yii\base\Controller-&gt;runAction()
#8 /var/www/html/prof-homepages/vendor/yiisoft/yii2/web/Application.php(103): yii\base\Module-&gt;runAction()
#9 /var/www/html/prof-homepages/vendor/yiisoft/yii2/base/Application.php(386): yii\web\Application-&gt;handleRequest()
#10 /var/www/html/prof-homepages/backend/web/index.php(16): yii\base\Application-&gt;run()
#11 {main}</pre>