R/empirical_estimate_PAFit.R

empirical_estimate_PAFit <- function(net, 
                                     outlier_prop = 95/100, deg_thresh = sqrt(max(final_deg))) {
  stats    <- get_statistics(as.PAFit_net(net), binning = FALSE,
                             only_PA = TRUE)
  result <- only_A_estimate(as.PAFit_net(net),stats, 
                            MLE = TRUE, stop_cond = 10^-5)
  

  A           <- result$estimate_result$A
  k           <- result$estimate_result$k
  names(A)    <- k
  names(k)    <- k
  var_new     <- result$estimate_result$var_A/(A[k == 1]^2)
  A           <- A/A[k == 1]
  log_var_new <- var_new/A^2
  upper_A     <- exp(log(A) + 2*sqrt(log_var_new))
  lower_A     <- exp(log(A) - 2*sqrt(log_var_new))
  names(upper_A) <- k
  names(lower_A) <- k
  names(var_new) <- k
  
  alpha_empi  <- result$estimate_result$alpha
  
  
  ci         <-  result$estimate_result$ci
  
  return(list(A = A, k = k, var = var_new, upper_A = upper_A, lower_A = lower_A,
              alpha = alpha_empi, ci = ci))
}
thongphamthe/mcPAFit documentation built on May 20, 2019, 10:23 p.m.