R/trans.armdata.r

In NMA: Network Meta-Analysis Based on Multivariate Meta-Analysis and Meta-Regression Models

trans.armdata <- function(study,treat1,treat2,n1,n2,y,SE,measure,data=NULL){

  data <- data.frame(data)

  if(measure=="logHR"){

	study <- data[, deparse(substitute(study))]
	treat1 <- data[, deparse(substitute(treat1))]
	treat2 <- data[, deparse(substitute(treat2))]
	n1 <- data[, deparse(substitute(n1))]
	n2 <- data[, deparse(substitute(n2))]
	y <- data[, deparse(substitute(y))]
	SE <- data[, deparse(substitute(SE))]

	study <- as.numeric(factor(study))

	trt <- c(treat1,treat2)

	treat <- levels(factor(trt))

	N <- max(study)
	p <- length(treat)
	
	L <- dim(data)[1]

	d1 <- d2 <- rep(NA,times=L)
	
	for(i in 1:L){
	
		yi <- y[i]
		vi <- SE[i]^2

		m1 <- n1[i]
		m2 <- n2[i]

		inve <- function(x){
		
			if((x[1]>0)&&(x[2]>0)){
	
				x1 <- x[1]/m1
				x2 <- x[2]/m2
			
				z1 <- log(-log(1-x1)) - log(-log(1-x2)) - yi
				z2 <- (((x1-1)*log(1-x1))^-2)*(x1*(1-x1)/m1) + (((x2-1)*log(1-x2))^-2)*(x2*(1-x2)/m2) - vi
			
				return(c(z1,z2))
			
			}

			if((x[1]<=0)||(x[2]<=0)){
	
				return(c(NaN,NaN))
			
			}
		
		}
		
		z0 <- c(.1*m1,.1*m2)
		d0 <- nleqslv(z0,inve)$x
		
		d1[i] <- d0[1]
		d2[i] <- d0[2]
	
	}
	
	id <- c(study,study)
	d <- c(d1,d2)
	n <- c(n1,n2)
	
	oi <- order(id,trt)

	Q1 <- data.frame(id,trt,d,n)[oi,]
	
	M <- dim(Q1)[1]

	j <- 1
	
	Q2 <- NULL
	
	repeat{
	
		wj <- which((Q1$id==Q1$id[j])&(Q1$trt==Q1$trt[j]))
	
		if(length(wj)==1){
			Q2 <- rbind(Q2,Q1[wj,])
			j <- j + 1
		}
		
		if(length(wj)!=1){
			qj <- Q1[wj[1],]
			qj$d <- mean(Q1$d[wj])
			qj$n <- mean(Q1$n[wj])
			Q2 <- rbind(Q2,qj)
			j <- max(wj) + 1
		}
		
		if(j>M)	break
	
	}

	colnames(Q2) <- c("study","trt","d","n")
	
	message("Note that the estimated numbers of events may not correspond to the actual event counts.")
	message("They are solely pseudo-data, designed so that the hazard ratio estimates generated by the setup function accord to the original data.")
	message("The NMA package tools rely solely on summary statistics for the synthesis analyses, so this is not problematic.")
	message("")

	return(Q2)

  }

  if(measure=="logOR"){

	study <- data[, deparse(substitute(study))]
	treat1 <- data[, deparse(substitute(treat1))]
	treat2 <- data[, deparse(substitute(treat2))]
	n1 <- data[, deparse(substitute(n1))]
	n2 <- data[, deparse(substitute(n2))]
	y <- data[, deparse(substitute(y))]
	SE <- data[, deparse(substitute(SE))]

	study <- as.numeric(factor(study))

	trt <- c(treat1,treat2)

	treat <- levels(factor(trt))

	N <- max(study)
	p <- length(treat)
	
	L <- dim(data)[1]

	d1 <- d2 <- rep(NA,times=L)
	
	for(i in 1:L){
	
		yi <- y[i]
		vi <- SE[i]^2

		m1 <- n1[i]
		m2 <- n2[i]

		inve <- function(x){
		
			if((x[1]>0)&&(x[2]>0)){
	
				x1 <- x[1]/m1
				x2 <- x[2]/m2
			
				z1 <- log( x[1]/(m1-x[1])/(x[2]/(m2-x[2])) ) - yi
				z2 <- 1/x[1] + 1/(m1-x[1]) + 1/x[2] + 1/(m2-x[2]) - vi
			
				return(c(z1,z2))
			
			}

			if((x[1]<=0)||(x[2]<=0)){
	
				return(c(NaN,NaN))
			
			}
		
		}
		
		z0 <- c(.1*m1,.1*m2)
		d0 <- nleqslv(z0,inve)$x
		
		d1[i] <- d0[1]
		d2[i] <- d0[2]
	
	}
	
	id <- c(study,study)
	d <- c(d1,d2)
	n <- c(n1,n2)
	
	oi <- order(id,trt)

	Q1 <- data.frame(id,trt,d,n)[oi,]
	
	M <- dim(Q1)[1]

	j <- 1
	
	Q2 <- NULL
	
	repeat{
	
		wj <- which((Q1$id==Q1$id[j])&(Q1$trt==Q1$trt[j]))
	
		if(length(wj)==1){
			Q2 <- rbind(Q2,Q1[wj,])
			j <- j + 1
		}
		
		if(length(wj)!=1){
			qj <- Q1[wj[1],]
			qj$d <- mean(Q1$d[wj])
			qj$n <- mean(Q1$n[wj])
			Q2 <- rbind(Q2,qj)
			j <- max(wj) + 1
		}
		
		if(j>M)	break
	
	}

	colnames(Q2) <- c("study","trt","d","n")
	
	message("Note that the estimated numbers of events may not correspond to the actual event counts.")
	message("They are solely pseudo-data, designed so that the odds ratio estimates generated by the setup function accord to the original data.")
	message("The NMA package tools rely solely on summary statistics for the synthesis analyses, so this is not problematic.")
	message("")

	return(Q2)

  }

  if(measure=="logRR"){

	study <- data[, deparse(substitute(study))]
	treat1 <- data[, deparse(substitute(treat1))]
	treat2 <- data[, deparse(substitute(treat2))]
	n1 <- data[, deparse(substitute(n1))]
	n2 <- data[, deparse(substitute(n2))]
	y <- data[, deparse(substitute(y))]
	SE <- data[, deparse(substitute(SE))]

	study <- as.numeric(factor(study))

	trt <- c(treat1,treat2)

	treat <- levels(factor(trt))

	N <- max(study)
	p <- length(treat)
	
	L <- dim(data)[1]

	d1 <- d2 <- rep(NA,times=L)
	
	for(i in 1:L){
	
		yi <- y[i]
		vi <- SE[i]^2

		m1 <- n1[i]
		m2 <- n2[i]

		inve <- function(x){
		
			if((x[1]>0)&&(x[2]>0)){
	
				x1 <- x[1]/m1
				x2 <- x[2]/m2
			
				z1 <- log( x1/x2 ) - yi
				z2 <- 1/x[1] - 1/m1 + 1/x[2] - 1/m2 - vi
			
				return(c(z1,z2))
			
			}

			if((x[1]<=0)||(x[2]<=0)){
	
				return(c(NaN,NaN))
			
			}
		
		}
		
		z0 <- c(.1*m1,.1*m2)
		d0 <- nleqslv(z0,inve)$x
		
		d1[i] <- d0[1]
		d2[i] <- d0[2]
	
	}
	
	id <- c(study,study)
	d <- c(d1,d2)
	n <- c(n1,n2)
	
	oi <- order(id,trt)

	Q1 <- data.frame(id,trt,d,n)[oi,]
	
	M <- dim(Q1)[1]

	j <- 1
	
	Q2 <- NULL
	
	repeat{
	
		wj <- which((Q1$id==Q1$id[j])&(Q1$trt==Q1$trt[j]))
	
		if(length(wj)==1){
			Q2 <- rbind(Q2,Q1[wj,])
			j <- j + 1
		}
		
		if(length(wj)!=1){
			qj <- Q1[wj[1],]
			qj$d <- mean(Q1$d[wj])
			qj$n <- mean(Q1$n[wj])
			Q2 <- rbind(Q2,qj)
			j <- max(wj) + 1
		}
		
		if(j>M)	break
	
	}

	colnames(Q2) <- c("study","trt","d","n")
	
	message("Note that the estimated numbers of events may not correspond to the actual event counts.")
	message("They are solely pseudo-data, designed so that the risk ratio estimates generated by the setup function accord to the original data.")
	message("The NMA package tools rely solely on summary statistics for the synthesis analyses, so this is not problematic.")
	message("")

	return(Q2)

  }


  if(measure=="RD"){

	study <- data[, deparse(substitute(study))]
	treat1 <- data[, deparse(substitute(treat1))]
	treat2 <- data[, deparse(substitute(treat2))]
	n1 <- data[, deparse(substitute(n1))]
	n2 <- data[, deparse(substitute(n2))]
	y <- data[, deparse(substitute(y))]
	SE <- data[, deparse(substitute(SE))]

	study <- as.numeric(factor(study))

	trt <- c(treat1,treat2)

	treat <- levels(factor(trt))

	N <- max(study)
	p <- length(treat)
	
	L <- dim(data)[1]

	d1 <- d2 <- rep(NA,times=L)
	
	for(i in 1:L){
	
		yi <- y[i]
		vi <- SE[i]^2

		m1 <- n1[i]
		m2 <- n2[i]

		inve <- function(x){
		
			if((x[1]>0)&&(x[2]>0)){
	
				x1 <- x[1]/m1
				x2 <- x[2]/m2
			
				z1 <- x1 - x2 - yi
				z2 <- x1*(1-x1)/m1 + x2*(1-x2)/m2 - vi
			
				return(c(z1,z2))
			
			}

			if((x[1]<=0)||(x[2]<=0)){
	
				return(c(NaN,NaN))
			
			}
		
		}
		
		z0 <- c(.1*m1,.1*m2)
		d0 <- nleqslv(z0,inve)$x
		
		d1[i] <- d0[1]
		d2[i] <- d0[2]
	
	}
	
	id <- c(study,study)
	d <- c(d1,d2)
	n <- c(n1,n2)
	
	oi <- order(id,trt)

	Q1 <- data.frame(id,trt,d,n)[oi,]
	
	M <- dim(Q1)[1]

	j <- 1
	
	Q2 <- NULL
	
	repeat{
	
		wj <- which((Q1$id==Q1$id[j])&(Q1$trt==Q1$trt[j]))
	
		if(length(wj)==1){
			Q2 <- rbind(Q2,Q1[wj,])
			j <- j + 1
		}
		
		if(length(wj)!=1){
			qj <- Q1[wj[1],]
			qj$d <- mean(Q1$d[wj])
			qj$n <- mean(Q1$n[wj])
			Q2 <- rbind(Q2,qj)
			j <- max(wj) + 1
		}
		
		if(j>M)	break
	
	}

	colnames(Q2) <- c("study","trt","d","n")
	
	message("Note that the estimated numbers of events may not correspond to the actual event counts.")
	message("They are solely pseudo-data, designed so that the risk difference estimates generated by the setup function accord to the original data.")
	message("The NMA package tools rely solely on summary statistics for the synthesis analyses, so this is not problematic.")
	message("")

	return(Q2)

  }


  if(measure=="MD"){

	study <- data[, deparse(substitute(study))]
	treat1 <- data[, deparse(substitute(treat1))]
	treat2 <- data[, deparse(substitute(treat2))]
	n1 <- data[, deparse(substitute(n1))]
	n2 <- data[, deparse(substitute(n2))]
	y <- data[, deparse(substitute(y))]
	SE <- data[, deparse(substitute(SE))]

	study <- as.numeric(factor(study))

	trt <- c(treat1,treat2)

	treat <- levels(factor(trt))

	N <- max(study)
	p <- length(treat)
	
	L <- dim(data)[1]

	d1 <- s1 <- rep(NA,times=L)
	
	for(i in 1:L){
	
		yi <- y[i]
		vi <- SE[i]^2

		m1 <- n1[i]
		m2 <- n2[i]

		inve <- function(x){
		
			if(x[2]>0){
	
				x1 <- x[1]
				x2 <- x[2]
			
				z1 <- (x1 - 0) - yi
				z2 <- x2*x2*((1/m1) + (1/m2)) - vi
			
				return(c(z1,z2))
			
			}

			if(x[2]<=0){
	
				return(c(NaN,NaN))
			
			}
		
		}
		
		z0 <- c(0,.1)
		d0 <- nleqslv(z0,inve)$x
		
		d1[i] <- d0[1]
		s1[i] <- d0[2]
	
	}
	
	mean0 <- rep(0,times=L)
	
	id <- c(study,study)
	m <- c(d1,mean0)
	s <- c(s1,s1)
	n <- c(n1,n2)
	
	oi <- order(id,trt)

	Q1 <- data.frame(id,trt,m,s,n)[oi,]
	
	M <- dim(Q1)[1]

	j <- 1
	
	Q2 <- NULL
	
	repeat{
	
		wj <- which((Q1$id==Q1$id[j])&(Q1$trt==Q1$trt[j]))
	
		if(length(wj)==1){
			Q2 <- rbind(Q2,Q1[wj,])
			j <- j + 1
		}
		
		if(length(wj)!=1){
			qj <- Q1[wj[1],]
			qj$m <- mean(Q1$m[wj])
			qj$s <- mean(Q1$s[wj])
			qj$n <- mean(Q1$n[wj])
			Q2 <- rbind(Q2,qj)
			j <- max(wj) + 1
		}
		
		if(j>M)	break
	
	}

	colnames(Q2) <- c("study","trt","m","s","n")
	
	message("Note that the estimated means and SDs may not correspond to the actual values.")
	message("They are solely pseudo-data, designed so that the mean difference estimates generated by the setup function accord to the original data.")
	message("The NMA package tools rely solely on summary statistics for the synthesis analyses, so this is not problematic.")
	message("")

	return(Q2)

  }

  # SMDは、seTEのほうが、n1, n2だけの関数になる。まあ、m1=d, m2=0, s=1 とすればよいのであるが、それであれば、逆算する必要はないよね。
  # そして、SMDが報告されていて、MDが報告されていない試験って、ないよね。
						
}