"calcmoment" <-
function(Y,time,deriv=0,k=1,central=F,realtime){
    N <- dim(Y)[2]
    moment <- NULL
    # Calculate the information function using the first deriv.
    I <- abs(t(apply(Y,2,diff)/diff(realtime)))
    time <- as.matrix(time[-1,])
    I <- I/apply(I,1,sum)
    I <- t(I/apply(I,1,sum))
    n <- dim(time)[1]
    # Compute the first k moments.
    for (j in 1:k)
      if ((j==1) | !central)
        moment <- rbind(moment,rep(1,n)%*%(I*(time^j)))
      else
        moment <- rbind(moment,rep(1,n)%*%(I*(t(t(time)-moment[1,])^j)))
    I <- rbind(I[1,],I)
    list(moments=moment,I=I)
  }

"criterion.mix" <-
function(theta,lambda.m,meanmoments,time,Y,q,b,target,lambda.t,smoothY,maxrange=0.8){
    K <- length(lambda.m)+1
    lambda.m <- lambda.m[c(1,1:(K-1))]
    lambda.x <- 1
    difftime <- diff(time)
    difftime <- difftime[c(1,1:length(difftime))]
    # alpha is the intercept or shift term.
    alpha <- theta[1]
    theta <- theta[-1]
    ef <-  as.vector(exp(b%*%theta)*difftime)
    # Compute X(t) less the intercept/shift term.
    int.ef <- cumsum(ef)+min(time)-ef[1]
    rawmoment <- calcmoment(as.matrix(Y),as.matrix(int.ef),deriv=1,k=1,realtime=time)$mom
    # Compute X(t).
    X <- alpha+int.ef
    # Calculate the moments of the synchronized curves.
    moment <- calcmoment(as.matrix(Y),as.matrix(X),deriv=1,k=K,central=T,realtime=time)
    moments <- moment$moments
    I <- as.vector(moment$I[-1])
      YofX <- predict(smooth.spline(X,Y),time)$y
    div <-(abs(max(X)-max(time))+abs(min(X)-min(time)))/(max(time)-min(time))
    penalty <- (10^3*(div-maxrange))^2*(div>maxrange)
    # Calculate the overall criterion including penalty for large
    # deviations from time end points.
    Q <-
      sum(lambda.m*((moments-meanmoments))^2)+lambda.x*sum((int.ef-time)^2) +
        lambda.t*(mean((YofX-target)^2)+penalty)
    Q}

"criterion.mom" <-
function(theta,lambda.m,meanmoments,time,Y,q,b,maxrange=0.8){
    K <- length(lambda.m)+1
    lambda.x <- 1
    # Compute X(t) less the intercept/shift term.
    difftime <- diff(time)
    difftime <- difftime[c(1,1:length(difftime))]
    ef <-  as.vector(exp(b%*%theta)*difftime)
    int.ef <- cumsum(ef)+min(time)-ef[1]
    rawmoment <-
  calcmoment(as.matrix(Y),as.matrix(int.ef),deriv=1,k=1,realtime=time)$mom
    # alpha is the intercept or shift term.
    alpha <- meanmoments[1]-rawmoment
    # Compute X(t).
    X <- alpha+int.ef
    # Calculate the moments of the synchronized curves.
    moment <-
      calcmoment(as.matrix(Y),as.matrix(X),deriv=1,k=K,central=T,realtime=time)
    moments <- moment$moments
    I <- as.vector(moment$I[-1])
    div <-(abs(max(X)-max(time))+abs(min(X)-min(time)))/(max(time)-min(time))
    penalty <- (10^3*(div-maxrange))^2*(div>maxrange)
    # Calculate the overall criterion including penalty for large
    # deviations from time end points.
    Q <- sum(lambda.m*((moments-meanmoments)[-1])^2+penalty)+lambda.x*sum((int.ef-time)^2)
    # Compute components of the derivative.
    b.ef <- apply(b*ef,2,cumsum)[-1,]
    b.ef2 <- apply(b*ef,2,cumsum)
    delalpha <- -matrix(apply(b.ef*I,2,sum),length(X)-1,q+1,byrow=T)
    delef1 <- matrix(ef[1]*b[1,],length(X),q+1,byrow=T)
    delX <- matrix(0,K-1,length(theta))
    for (k in 2:K)
      delX[k-1,]  <-k*apply((X[-1]-meanmoments[1])^(k-1)*(b.ef+delalpha)*I,2,sum)
    # Compute the derivative.
    delQ <-
      2*as.vector((lambda.m*((moments-meanmoments)[-1])))%*%delX+2*lambda.x*apply((b.ef2-delef1)*(int.ef-time),2,sum)
    attr(Q,"gradient") <- delQ
    Q}

"curvesync" <-
function(Y,realtime=seq(0,1,length=dim(Y)[1]),k=3,q=5,trans=matrix(0,N,q+2),plotgraphs=T,lambda=100,diags=F,linear=T,df=10,shift=F,CAM=T,lambda.p=NULL,target=NULL,maxit=100,maxrange=0.8,cs.target=F){
    library(splines)
    options(show.error.messages=F)
    N <- dim(Y)[2]
    if (length(lambda)==1)
      lambda <- rep(lambda,k-1)
    cam <- camp <- NULL
    # Call CAM function if CAM=T.
    if (CAM)
      cam <-
        curvesync.mom(Y=Y,realtime=realtime,k=k,q=q,trans=trans,plotgraphs=plotgraphs,lambda=lambda,diags=diags,linear=linear,df=df,shift=shift,maxrange=maxrange)
    # Call CAMP function if lambda.p is not NULL.
    if (!is.null(lambda.p))
      camp <- curvesync.mix(Y=Y,realtime=realtime,k=k,q=q,trans=trans,plotgraphs=plotgraphs,lambda=lambda,diags=diags,df=df,lambda.t=lambda.p,target=target,maxit=maxit,maxrange=maxrange,cs.target=cs.target)
    list(cam=cam,camp=camp)}

"curvesync.mix" <-
function(Y,realtime=seq(0,1,length=dim(Y)[1]),k=3,q=5,trans=matrix(0,N,q+2),plotgraphs=F,lambda=100,diags=F,df=10,lambda.t=10,target=NULL,maxit=100,maxrange=maxrange,cs.target=F){
    N <- dim(Y)[2]
    varx <- 0
    smoothY <- as.list(rep(0,N))
    # Smooth the curves.
    for (i in 1:N){
      smoothY[[i]] <- smooth.spline(realtime,Y[,i],df=df)
      Y[,i] <- smoothY[[i]]$y}
    # Compute target function.
    if (is.null(target))
      if (cs.target)
        target <- smooth.spline(apply(Y,1,mean))$y
      else
        target <-apply(curvesync.mom(Y,realtime=realtime,k=k,q=q,plotgraphs=F,lambda=100,diags=F,linear=T,df=df,shift=F)$YofX,1,mean)
    n <- length(realtime)
    meancurve <- varcurve <- rep(0,n)
    meanmoments <- rep(0,k)
    time <- matrix(realtime,n,N)
    par(mfrow=c(2,2))
    h <- diff(range(realtime))/n
    YofX <- Y
    origmeancurve <- apply(Y,1,mean)
    # Plot the smoothed but unsynchronized curves.
    if (plotgraphs){
      plot(realtime,realtime,ylim=range(Y),type='n',
          ylab="Observed curves", xlab="Actual time")
      for (i in 1:N)
          lines(realtime,Y[,i])
      lines(realtime,origmeancurve,col=2,lwd=4)
      plot(realtime,realtime,ylim=range(Y),type='n',
    ylab="Synchronized curves", xlab="Synchronized time")
      }
    time <- matrix(realtime,n,N)
    # Compute the first k moments on the smoothed data.
    moments <- calcmoment(Y,time,k=k,central=T,realtime=realtime)
    I <- moments$I
    mom <- moments$moments
    # Compute the estimated first k moments for the generating function Z.
    for (i in 1:k)
      meanmoments[i] <- mean(sign(mom[i,])*(abs(mom[i,])^(1/i)))^i
    moments <- mom
    # Standardize lambda.
    if (k==2)
      lambda <- lambda/mean(moments[2,])
    else
      lambda <- lambda/apply(abs(moments[-1,]),1,mean)
    # Create basis for synchronization functions.
    b <- bs(realtime,df=q+1,intercept=T)
    # This loops calculates the synchronization functions.
    for (i in 1:N){
      smax <- 10
      tempfit <- NULL
      # Call CAMP optimizer.
      while (is.null(tempfit$est)){
        tempfit <- try(nlm(f=criterion.mix, p=trans[i,],
                           check.analyticals=T, lambda=lambda,
                           meanmoments=meanmoments,time=realtime,Y=Y[,i],q=q,b=b,target=target,lambda.t=lambda.t,smoothY=smoothY[[i]],maxrange=maxrange,print.level=0,iterlim=maxit,ndigit=8,stepmax=smax))
        smax <- smax/2}
      # Synchronization function coefficients.
      trans[i,] <- tempfit$est
      # Compute X(t) minus the shift term.
      difftime <- diff(realtime)
      difftime <- difftime[c(1,1:length(difftime))]
      ef <-  as.vector(exp(b%*%trans[i,-1])*difftime)
      int.ef <- cumsum(ef)+min(realtime)-ef[1]
      if (diags)
        print(paste("Optimization code for function ",i," = ",tempfit$code))
      # Compute X(t) including the shift term
      time[,i] <- trans[i,1]+int.ef
      # Plot the synchronized curves.
      if (plotgraphs)
      lines(time[,i],Y[,i])
    }
    # Compute the synchronizatied curves. 
    for (i in 1:N)
      YofX[,i] <- predict(smooth.spline(time[,i],Y[,i],df=df),realtime)$y
    meancurve <- apply(YofX,1,mean)
    # Plot the synchronization functions.
    if (plotgraphs){
      lines(realtime,meancurve,col=2,lwd=4)
     plot(realtime,realtime,xlim=range(time),type='n',ylab=
               "Actual time",xlab="Warped time") 
    for (i in 1:N)
        lines(time[,i],realtime)
    abline(0,1,lwd=4,col=2)
    }
    # Compute quality of synchronization statistic.
    sync <- 0
    for (i in 1:N){
      origmeancurve <- apply(Y[,-i],1,mean)
      meancurve <- apply(YofX[,-i],1,mean)
      sync <- sync+ sum((YofX[,i]-meancurve)^2)/sum((Y[,i]-origmeancurve)^2)}
    list(trans=trans,YofX=YofX,X=time,mom=mom,Y=Y,sync=sync/N)}

"curvesync.mom" <-
function(Y,realtime=seq(0,1,length=dim(Y)[1]),k=3,q=5,trans=matrix(0,N,q+2),plotgraphs=F,lambda=100,diags=F,linear=T,df=10,shift=F,maxrange=0.8){
    N <- dim(Y)[2]
    smoothY <- as.list(rep(0,N))
    # Smooth the curves.
    for (i in 1:N){
      smoothY[[i]] <- smooth.spline(realtime,Y[,i],df=df)
      Y[,i] <- smoothY[[i]]$y}
    n <- length(realtime)
    meancurve <- varcurve <- rep(0,n)
    meanmoments <- rep(0,k)
    time <- matrix(realtime,n,N)
    par(mfrow=c(2,2))
    h <- diff(range(realtime))/n
    YofX <- Y
    origmeancurve <- apply(Y,1,mean)
    # Plot the smoothed but unsynchronized curves.
    if (plotgraphs){
      plot(realtime,realtime,ylim=range(Y),type='n',
   ylab="Observed curves", xlab="Actual time")
      for (i in 1:N)
          lines(realtime,Y[,i])
      lines(realtime,origmeancurve,col=2,lwd=4)
      plot(realtime,realtime,ylim=range(Y),type='n',
    ylab="Synchronized curves", xlab="Synchronized time")
    }
    time <- matrix(realtime,n,N)
    # Compute the first k moments on the smoothed data.
    moments <- calcmoment(Y,time,k=k,central=T,realtime=realtime)
    I <- moments$I
    mom <- moments$moments
    # Compute the estimated first k moments for the generating function Z.
    for (i in 1:k)
      meanmoments[i] <- mean(sign(mom[i,])*(abs(mom[i,])^(1/i)))^i
    moments <- mom
    # Standardize lambda.
    if (k==2)
      lambda <- lambda/mean(moments[2,])
    else
      lambda <- lambda/apply(abs(moments[-1,]),1,mean)
    # Create basis for synchronization functions.
    b <- bs(realtime,df=q+1,intercept=T)
    # This loops calculates the synchronization functions.
    for (i in 1:N){
      # Initialize the synchronization function coefficients.
      if (trans[i,1]==0){
        beta <- sqrt(meanmoments[2]/moments[2,i])
        if (shift)
          beta <- 1
        trans[i,-1] <-rep(log(beta),q+1)
        # Call non-linear optimizer if doing non-linear CAM.
        if (!linear & !shift){
          tempfit <- nlm(f=criterion.mom, p=trans[i,-1],
                         check.analyticals=F, lambda=lambda,
                         meanmoments=meanmoments,time=realtime,Y=Y[,i],q=q,b=b,maxrange=maxrange,print.level=0,iterlim=200)
          # Synchronization function coefficients.
          trans[i,-1] <- tempfit$est
        }
        # Compute X(t) minus the shift term.
        difftime <- diff(realtime)
        difftime <- difftime[c(1,1:length(difftime))]
        ef <-  as.vector(exp(b%*%trans[i,-1])*difftime)
        int.ef <- cumsum(ef)+min(realtime)-ef[1]
        # Compute the shift term.
        trans[i,1] <-mean(moments[1,])-calcmoment(as.matrix(Y[,i]),as.matrix(int.ef),k=1,central=T,realtime=realtime)$mom
        if (diags){
          print(paste("Optimization code for function ",i," = ",tempfit$code))
        }}
      # Compute X(t) including the shift term
      time[,i] <- trans[i,1]+int.ef
      # Plot the synchronized curves.
      if (plotgraphs)
      lines(time[,i],Y[,i])
    }
    # Compute the synchronizatied curves.
    for (i in 1:N)
      YofX[,i] <- predict(smooth.spline(time[,i],Y[,i],df=df),realtime)$y
    meancurve <- apply(YofX,1,mean)
    # Plot the synchronization functions.
    if (plotgraphs){
    lines(realtime,meancurve,col=2,lwd=4)                    
    plot(realtime,realtime,xlim=range(time),type='n',ylab=
               "Actual time",xlab="Warped time")
    for (i in 1:N)
        lines(time[,i],realtime)
    abline(0,1,lwd=4,col=2)
  }
    # Compute quality of synchronization statistic.
    sync <- 0
    for (i in 1:N){
      origmeancurve <- apply(Y[,-i],1,mean)
      meancurve <- apply(YofX[,-i],1,mean)
      sync <- sync+ sum((YofX[,i]-meancurve)^2)/sum((Y[,i]-origmeancurve)^2)}
    list(trans=trans,YofX=YofX,X=time,mom=mom,Y=Y,sync=sync/N)}

"testY" <-
structure(c(8.48758852272934e-08, 1.04767104491598e-06, 1.03335596846737e-05, 
8.14440477419407e-05, 0.000512924134210825, 0.00258125734665798, 
0.0103799165924274, 0.0333534540868078, 0.0856391212487582, 0.175706686364499, 
0.288065232655057, 0.377391091775398, 0.396234222946465, 0.348781049889672, 
0.262320948664376, 0.171413007517112, 0.104865608867808, 0.0760200826734867, 
0.0840195616155796, 0.118292475637086, 0.161429330617292, 0.193057451324636, 
0.197631453535480, 0.172336610286234, 0.127886560401772, 0.0807444509769802, 
0.0433734294057039, 0.0198223301489192, 0.00770734552360848, 
0.00254961142000117, 0.000717566537195065, 0.000171818280934472, 
3.50022183793077e-05, 6.06653599843951e-06, 8.94550775233508e-07, 
1.12224640337849e-07, 1.19781695573824e-08, 1.08770643808051e-09, 
8.40334028387313e-11, 5.52345995395264e-12, 3.08879778753422e-13, 
1.46955807618188e-14, 5.94843935912243e-16, 2.04851114788092e-17, 
6.00195359450288e-19, 1.49611832793931e-20, 3.17291371246328e-22, 
5.72492219052244e-24, 8.78819983141138e-26, 1.14775519521532e-27, 
9.96567785067151e-08, 1.07767663547691e-06, 9.51630438263597e-06, 
6.86193210014806e-05, 0.000404038492712436, 0.00194266225711328, 
0.00762728716642305, 0.0244535058195823, 0.0640192265907753, 
0.136860491504847, 0.238915807835518, 0.340576957354821, 0.396475411760247, 
0.283420372172378, 0.0845325194263037, 0.130984531074709, 0.196980669775694, 
0.0890383605309067, 0.0117589468067520, 0.000453575538415081, 
5.10997768650325e-06, 1.68142128948288e-08, 1.61592886706415e-11, 
4.53582296617501e-15, 3.71858840726895e-19, 8.90406425770327e-24, 
6.2271111690648e-29, 1.27195899693863e-34, 7.58835715287115e-41, 
1.32224086053062e-47, 6.72917067289878e-55, 1.00023183332843e-62, 
4.34238134326122e-71, 5.5060912808428e-80, 2.03913985409225e-89, 
2.20566106912578e-99, 6.96816942737654e-110, 6.42964189003555e-121, 
1.73277840524453e-132, 1.36391538438902e-144, 3.13559518059922e-157, 
2.10542951843100e-170, 4.12904471282334e-184, 2.36508428917547e-198, 
3.95668812502057e-213, 1.93332637164991e-228, 2.75909522436711e-244, 
1.15004904862700e-260, 1.40008433229054e-277, 4.97828994665947e-295, 
3.46378699302355e-22, 2.3783037013861e-19, 1.01832308785226e-16, 
2.71897999905276e-14, 4.52718851668301e-12, 4.70059992866691e-10, 
3.04354791906737e-08, 1.22887905452653e-06, 3.09414230685082e-05, 
0.000485817964462101, 0.00475673547054375, 0.0290433414865476, 
0.110582501845507, 0.262560886139198, 0.388777819757844, 0.381260901233215, 
0.300843033943979, 0.196201449974100, 0.113266505823635, 0.0763789229966975, 
0.0862668198738704, 0.127617584993475, 0.174591323622528, 0.199107512432707, 
0.185670054969778, 0.141119948310887, 0.0873772154921127, 0.0440690603477738, 
0.0181046051810237, 0.00605847059295754, 0.00165141210924259, 
0.000366662430873419, 6.63125678409547e-05, 9.76885718766093e-06, 
1.17222278449466e-06, 1.14576364804468e-07, 9.12217506433194e-09, 
5.91589216048566e-10, 3.12507595743093e-11, 1.34468137444309e-12, 
4.71299179686718e-14, 1.34552713675482e-15, 3.12900834278379e-17, 
5.92706191353677e-19, 9.14514697125126e-21, 1.14937143012974e-22, 
1.17665349400549e-24, 9.81194855719449e-27, 6.66469714584513e-29, 
3.68743232534399e-31, 2.10185100752548e-12, 4.42424474387962e-11, 
7.70076158770194e-10, 1.10837196050341e-08, 1.3191515771783e-07, 
1.29825854464948e-06, 1.05653769484118e-05, 7.10994303247426e-05, 
0.000395643981629497, 0.00182054041618638, 0.00692713105757895, 
0.0217953478750751, 0.0567063301992184, 0.121999121335502, 0.217040000598027, 
0.31928874870018, 0.388422531900997, 0.390632913095709, 0.32402688968585, 
0.223182197882086, 0.133002944457514, 0.0825563619246442, 0.0783639220257563, 
0.110375877972829, 0.157517181742928, 0.193119129738817, 0.196344665143482, 
0.164501452198470, 0.113450126413242, 0.0643937857308306, 0.0300796525680676, 
0.0115634969006046, 0.00365841248838133, 0.0009525399385604, 
0.000204108337787155, 3.59936031922172e-05, 5.22368439058969e-06, 
6.23900826117368e-07, 6.13255170037104e-08, 4.96082494337256e-09, 
3.30258201937115e-10, 1.80942400834653e-11, 8.15857125862409e-13, 
3.02743550581286e-14, 9.24533634171446e-16, 2.32357960980236e-17, 
4.80595195122883e-19, 8.1806595374012e-21, 1.14599900065761e-22, 
1.32119580338083e-24, 5.86230012509651e-07, 1.20244426306749e-05, 
0.000167932838690903, 0.00159690864721246, 0.0103394701996313, 
0.0455816752303972, 0.136822091067902, 0.279639187819263, 0.389169385363716, 
0.384708626744132, 0.314133368344105, 0.216417090876316, 0.131059734609074, 
0.0829062049257119, 0.0775409659537045, 0.106402989876032, 0.151377439247225, 
0.188911490161273, 0.198823941221753, 0.175144173637697, 0.128952198869462, 
0.0793333135427032, 0.0407807880907633, 0.0175156005094959, 0.00628583752585637, 
0.00188481984961176, 0.000472220780261417, 9.88527806043893e-05, 
1.72902483991437e-05, 2.5268654779691e-06, 3.08554268365642e-07, 
3.14810741250905e-08, 2.68371301414753e-09, 1.91157415185274e-10, 
1.13766598537881e-11, 5.65727049525338e-13, 2.35054008794461e-14, 
8.16012717390757e-16, 2.36698082069472e-17, 5.73668624186488e-19, 
1.16170538655357e-20, 1.96561951629235e-22, 2.77889063308355e-24, 
3.28255492146362e-26, 3.2398215833315e-28, 2.67176794108089e-30, 
1.84096241207860e-32, 1.05988721939041e-34, 5.09850601572796e-37, 
2.04924915071852e-39, 1.12553969335275e-25, 6.56142557217746e-23, 
2.61461537969379e-20, 7.12179763649254e-18, 1.32600041518636e-15, 
1.6876015052278e-13, 1.46814300834558e-11, 8.7304982889512e-10, 
3.54880319141572e-08, 9.86045627584018e-07, 1.87276947241600e-05, 
0.00024313296388885, 0.00215762308730029, 0.0130881882133233, 
0.0542694870460104, 0.153817007724432, 0.298008113287559, 0.394691629462382, 
0.341054562544846, 0.176228567411742, 0.0789953025314174, 0.103091522809636, 
0.178443680588552, 0.193025551491316, 0.121454788545845, 0.0442793293581995, 
0.00935158634047718, 0.00114409683821673, 8.10836493170355e-05, 
3.32886957027268e-06, 7.91686542104526e-08, 1.09069352826646e-09, 
8.7045274785802e-12, 4.02421018229628e-14, 1.07772783862532e-16, 
1.67197791810074e-19, 1.50260546549229e-22, 7.82262265817359e-26, 
2.35913497691367e-29, 4.12141348898448e-33, 4.17092688969455e-37, 
2.44518631865436e-41, 8.30394239043873e-46, 1.63361490990355e-50, 
1.86169288998595e-55, 1.22902122267810e-60, 4.70006338775808e-66, 
1.04121649238363e-71, 1.33619964645294e-77, 9.93332527640905e-84, 
6.23774477036953e-23, 8.97548366679356e-21, 9.9693652944993e-19, 
8.54784319447037e-17, 5.65750533569411e-15, 2.89049875215666e-13, 
1.13998723587633e-11, 3.47061611989768e-10, 8.15628812385577e-09, 
1.47964634887866e-07, 2.07206169914565e-06, 2.23989070670222e-05, 
0.000186909063333722, 0.00120396339470371, 0.00598653593204763, 
0.0229782683286930, 0.068082958256233, 0.155718121383579, 0.274929300209186, 
0.374710306771120, 0.391780061951666, 0.285230513393455, 0.145859919314482, 
0.077117575810229, 0.0979675877880103, 0.163390095277457, 0.199595042172788, 
0.164668286282432, 0.0911211145016391, 0.0338016324805844, 0.00840517679838213, 
0.00140102226206695, 0.000156542441859238, 1.17248702821091e-05, 
5.88671300775357e-07, 1.98119194775861e-08, 4.4696047416267e-10, 
6.75928187950775e-12, 6.8520554848072e-14, 4.65618318458801e-16, 
2.12093871278673e-18, 6.4761224891255e-21, 1.32553397479936e-23, 
1.81867688450130e-26, 1.67266551901120e-29, 1.03122049470045e-32, 
4.26169967256205e-36, 1.18060065735799e-39, 2.19236103283825e-43, 
2.7290406327682e-47, 6.21384000382822e-21, 4.59839808838576e-18, 
2.02987635995522e-15, 5.34501743655165e-13, 8.3954733521095e-11, 
7.86606931965311e-09, 4.39629545598612e-07, 1.46565750072975e-05, 
0.000291470521974980, 0.00345758848712091, 0.0244663111146225, 
0.103271454895250, 0.260022074855611, 0.390557571252522, 0.373967186480595, 
0.275007663887108, 0.161162163390228, 0.0889490730236073, 0.0775348337327268, 
0.114660344689396, 0.168334937769442, 0.198896463336412, 0.182680308201123, 
0.129812762070773, 0.0713235886761868, 0.0302972700624519, 0.00995001283626282, 
0.0025263469881011, 0.000495919851393606, 7.52625815718804e-05, 
8.83071464898794e-06, 8.01053854629813e-07, 5.61793487349006e-08, 
3.04607651737774e-09, 1.27689162801957e-10, 4.13824684868695e-12, 
1.03687849940664e-13, 2.0085769187832e-15, 3.00814100022410e-17, 
3.48302853833603e-19, 3.11792022625181e-21, 2.15785415552564e-23, 
1.15459135123876e-25, 4.77621242639609e-28, 1.52752417066442e-30, 
3.77695378131455e-33, 7.22011927959832e-36, 1.06707807570588e-38, 
1.21926229613023e-41, 1.07707814535163e-44, 2.49520693373523e-32, 
1.26231576305125e-28, 3.72630519273413e-25, 6.41856313349874e-22, 
6.45128446576564e-19, 3.78358544937616e-16, 1.29482143512651e-13, 
2.58562294683825e-11, 3.01279433319748e-09, 2.04843768811829e-07, 
8.12690173780448e-06, 0.000188137877106131, 0.00254142003518663, 
0.0200320613900772, 0.0921348396365573, 0.247270559369179, 0.387252753998905, 
0.371056851743868, 0.254504610117424, 0.133259754298284, 0.0766981029360448, 
0.0941939082826497, 0.152608745032958, 0.196657860902473, 0.184376395815914, 
0.124550929899811, 0.060563047339623, 0.0211955305720029, 0.00533890667105733, 
0.000967902004016371, 0.000126293624715086, 1.18604848024262e-05, 
8.01667779306625e-07, 3.89993519189302e-08, 1.36549880808656e-09, 
3.44109160934228e-11, 6.24125774276205e-13, 8.14739880479305e-15, 
7.65485213561844e-17, 5.17637357091756e-19, 2.5193319254676e-21, 
8.8250290391932e-24, 2.22493777758252e-26, 4.03729580875649e-29, 
5.27271437094183e-32, 4.95620081192671e-35, 3.35300744109784e-38, 
1.63264256125737e-41, 5.72161570078352e-45, 1.44316912522926e-48, 
1.47269629797484e-15, 3.37343983735427e-14, 6.61837196231151e-13, 
1.11211014321171e-11, 1.60052646709817e-10, 1.97286035991372e-09, 
2.08280332031526e-08, 1.8832960868189e-07, 1.45850316123875e-06, 
9.67418913098879e-06, 5.49591744571391e-05, 0.000267414071010338, 
0.00111441446037816, 0.00397765992890681, 0.0121598156553126, 
0.0318379257476735, 0.0713971696373322, 0.137131045202565, 0.225584701799820, 
0.317837413320259, 0.383560364526581, 0.392358696756109, 0.284799011106194, 
0.143764323883646, 0.0764756752998113, 0.100704508980012, 0.167405124830604, 
0.199444940284275, 0.158319057015887, 0.0832407232244305, 0.0289756847042632, 
0.00667746563378314, 0.00101875217550312, 0.000102897371595973, 
6.88048181248543e-06, 3.04587627638034e-07, 8.92656768993475e-09, 
1.73195093277857e-10, 2.22466950993551e-12, 1.89179580655736e-14, 
1.06503002496230e-16, 3.9694312239606e-19, 9.79430330491956e-22, 
1.59991773806727e-24, 1.73021708071550e-27, 1.23874631026311e-30, 
5.87141786480044e-34, 1.84239371359678e-37, 3.82737251207934e-41, 
5.2637815924726e-45), .Dim = as.integer(c(50, 10)))