metaanalysis.R

## meta analysis, sample size based
metaZn =  function(zdisc,zrepl,ndisc,nrepl)
  {
  ## calculate meta analysis Zscore
  ## zdisc and zrepl are zscores in the discovery and replication sets respectively
  ## ndisc and nrepl are sample sizes in the discovery and replication sets
    wdisc = sqrt(ndisc)
    wrepl = sqrt(nrepl)
    ( zdisc * wdisc  + zrepl * wrepl )/sqrt( wdisc^2 + wrepl^2 )
}

metaZn3 =  function(zdisc,zrepl,zrepl2,,ndisc,nrepl,nrepl2)
  {
  ## calculate meta analysis Zscore
  ## zdisc and zrepl are zscores in the discovery and replication sets respectively
  ## ndisc and nrepl are sample sizes in the discovery and replication sets
    wdisc = sqrt(ndisc)
    wrepl = sqrt(nrepl)
    wrepl2 = sqrt(nrepl2)
    ( zdisc * wdisc  + zrepl * wrepl + zrepl2 * wrepl2)/sqrt( wdisc^2 + wrepl^2 + wrepl2^2)
}

## this is robust to correlation between entries
## it uses the fact that the average of Cauchy r.v. is Cauchy regardless of correlation between them
## https://www.cell.com/ajhg/pdfExtended/S0002-9297(19)30002-3
acat = function(pvec) 
{
  TT = sum( tan( (0.5 - pvec) *pi ) )
  .5 - atan(TT / length(pvec)) / pi
}
## find better implementation here https://github.com/yaowuliu/ACAT/blob/7f7d43162b098ad72315793bc1020e00e48c043d/R/ACAT.R#L80
## if zscores have different signs use sacat instead of acat so that opposite signs get cancelled out rather than increase significance
sacat = function(pvec,svec=1) 
{
  if(abs(sum(svec)) != sum(abs(svec)) )
    {
      print("some of the entries have opposite signs")
      TTplus = sum( svec * tan( (0.5 - pvec) *pi ) )
      TTminus = - sum( svec * tan( (0.5 - pvec) *pi ) )
      pplus = .5 - atan(TTplus / length(pvec)) / pi
      pminus = .5 - atan(TTminus / length(pvec)) / pi
      print(pplus)
      print(pminus)
      pacat = 2 * min(pplus, pminus)
      sacat = ifelse( pacat == pplus *2, 1, -1)
      print(pacat)
      print(sacat) 
      res = pacat
      } else
    {
      print("all signs are the same, using acat(pvec)")
      res = acat(pvec) 
    }
    res
  ## output sign_acat * pacat
}



## fisher's meta analysis, does not take into account direction of effect
## assumes independents of entries
## https://en.wikipedia.org/wiki/Fisher%27s_method
fisher = function(pvec)  
{
  X2 = -2 * sum(log(pvec))
  1 - pchisq(X2, 2*length(pvec))
}