Curvas reales, media y simuladas
```r library(iAR) library(polynom) library(dplyr) library(readr) library(simex) library(xtable) require(Rdrw)
setwd(“C:/Users/cataluna/OneDrive/Escritorio/usach/Tesis/AGN/AGN”)
AGN1=read.table(“AGN1.csv”,sep=“,”,header=T) AGN2=read.table(“AGN2.csv”,sep=“,”,header=T) AGN3=read.table(“AGN3.csv”,sep=“,”,header=T) AGN4=read.table(“AGN4.csv”,sep=“,”,header=T) AGN5=read.table(“AGN5.csv”,sep=“,”,header=T) AGN6=read.table(“AGN6.csv”,sep=“,”,header=T) AGN7=read.table(“AGN7.csv”,sep=“,”,header=T) AGNF=rbind(AGN1,AGN2,AGN3,AGN4,AGN5,AGN6,AGN7)
aggregate(mjd~fid,AGNF,range)
object=names(table(AGNF[,1])) l0=length(object) crossmatchedobjects=read.csv(“ts_v9.0.1_SMBH_ZTF_xmatch.csv”) object=data.frame(“oid”=object) join=merge(object, crossmatchedobjects, by = “oid”)
MSEMLEg=rep(0,l0) o2=iAR::utilities() o2<-gentime(o2, n=400, distribution = “expmixture”, lambda1 = 55, lambda2 = 3,p1 = 0.15, p2 = 0.85)
newdata=NULL simdata=NULL
for(i in 1:l0) { obj=join\(oid[i] class=join\)survey_class_mapped[i] pos=which(AGNF[,1]==obj) lc=AGNF[pos,] lcg=na.omit(lc[which(lc[,3]==1),c(2,4,5)]) lcr=na.omit(lc[which(lc[,3]==2),c(2,4,5)]) lcg=lcg[order(lcg[,1]),] pos=which(diff(lcg[,1])==0)
if(length(pos)>0) { for(k in pos) { lcg[k,]=apply(lcg[c(k,k+1),],2,mean) } lcg=lcg[-c(pos+1),] }
p=which(lcg[,2]>=30) if(length(p)>0) lcg=lcg[-p,]
if(dim(lcg)[1]<10) next
p=which(lcr[,2]>=30) if(length(p)>0) lcr=lcr[-p,]
lcr=lcr[order(lcr[,1]),] pos=which(diff(lcr[,1])==0)
if(length(pos)>0) { for(k in pos) { lcr[k,]=apply(lcr[c(k,k+1),],2,mean) } lcr=lcr[-c(pos+1),] }
if(dim(lcr)[1]<10) next
media=mean(lcg[,2]) primervalor=c(58242.22,media,0.01) ultimovalor=c(60316.16,media,0.01) lcg=rbind(primervalor,lcg,ultimovalor)
lcg <- lcg[order(lcg[,1]), , drop = FALSE] lcg <- lcg[!duplicated(lcg[,1]), , drop = FALSE] d <- diff(lcg[,1])
if (any(d <= 0)) { eps <- .Machine$double.eps^0.5 for (j in which(d <= 0)) lcg[j+1,1] <- lcg[j,1] + eps }
fid=1 new=data.frame(obj,lcg[,1],fid,lcg[,2],lcg[,3]) names(new)=names(AGNF) newdata=rbind(newdata,new)
n=dim(lcg)[1] x <- iAR(family = “norm”, times = lcg[,1], series = lcg[,2], series_esd = lcg[,3]) x = loglik(x)
st <- sort(unique(as.numeric(o2@times))) simiar <- iAR(family = “norm”, times = st, coef = x@coef) simiar <- iAR::sim(simiar) newsim <- data.frame(obj, st, fid, simiar@series + media, 0)
names(newsim)=names(AGNF) simdata=rbind(simdata,newsim) }