martinig / Fitness-and-dispersal-MA

https://martinig.github.io/Fitness-and-dispersal-MA/
0 stars 0 forks source link

how to convert percent/proportion alive #3

Closed martinig closed 6 months ago

martinig commented 10 months ago
itchyshin commented 10 months ago

There are a couple of parts

  1. first we need to proportion into logit scale (also need to know what to do with when 0 and 100%)
  2. Then, we turn into SMD
  3. Then to Zr via r

I think 2 and 3 are covered by one of the functions I sent so you can see whether you can find out

for the first step use logit() in the car package (which resolves the issue of 0 and 100%)

smd = (logit(p1) - logit(p2))/ (pi/sqrt(3)) # this is step 2 (p1 is proportion for 1 group and p2 for the other)

martinig commented 9 months ago

Here are three examples.

Example 1: N1 = 1 N2 = 13 prop1 =1.00 prop2=0.9230769231

Example 2: N1 = 256 N2 = 157 prop1 = 0 prop2=0.4840764331

Example 3: N1 = 8 N2 = 0 prop1 =0.25 prop2=NA

martinig commented 9 months ago
#function to convert proportion to r value 
# calculate_SMD takes two proportions (p1 and p2) as input, converts them into logit scale using the logit function from the car package, and then calculates the SMD using the provided formula.

library(car) # Load the car package for the logit function 

calculate_smd <- function(p1, p2, n1, n2) {

 # Step 1: Convert proportions to logit scale
  logit_p1 <- logit(p1)
  logit_p2 <- logit(p2)

  # Step 2: Calculate the Standardized Mean Difference (SMD)
  smd <- (logit_p1 - logit_p2) / (pi / sqrt(3))

  # Step 3: Adjust SMD for sample sizes
  smd_adjusted <- smd * sqrt((n1 + n2) / (n1 * n2))

  return(smd_adjusted)
}

#example 1:
n1=1
n2=13
p1=1
p2=0.9230769231

result <- calculate_smd(p1, p2, n1, n2) #gives me a warning message
print(result) #0.6743569

#example 2:
n1=256
n2=157
p1=0
p2=0.4840764331

result <- calculate_smd(p1, p2, n1, n2) #also gives me a warning message
print(result) #-0.201187

What to do for example 3?

I'm also not sure how to do step 3 in your original response.

itchyshin commented 9 months ago

@martinig - how to deal with 0 and 1 is case-by-case - we can leave this till we need to analyses data (you will find out why I leave this later)

itchyshin commented 6 months ago

sorted - we use arcsine transformation

martinig commented 6 months ago

This is the code for this:

}else if(method == "proportion_method1"){

    h <- n/n1 + n/n2
    p <- n1/n # prop for n1
    q <- n2/n # prop for n2
    s1 <- 1/sqrt(8)
    s2 <- 1/sqrt(8)
    m1 <- asin(sqrt(m1))
    m2 <- asin(sqrt(m2))
    s_pool <- sqrt( ((n1-1)*s1^2 + (n2-1)*s2^2) / (n - 2) )
    j <- 1 - (3 / (4*n - 9))
    d <- ((m2 - m1) / s_pool) * j
    r_pb <-  d / sqrt(d^2 + h)
    r_b <- r_pb*(sqrt(p*q)/dnorm(qnorm(p)))
    r <- r_b #r_b = r 

  }else if(method == "proportion_method2"){

    h <- n/n1 + n/n2
    p <- n1/n # prop for n1
    q <- n2/n # prop for n2
    s1 <- 1/sqrt(sd1^2/(4*m1*(1-m1)))
    s2 <- 1/sqrt(sd2^2/(4*m2*(1-m2)))
    m1 <- asin(sqrt(m1/100))
    m2 <- asin(sqrt(m2/100))
    s_pool <- sqrt( ((n1-1)*s1^2 + (n2-1)*s2^2) / (n - 2) )
    j <- 1 - (3 / (4*n - 9))
    d <- ((m2 - m1) / s_pool) * j
    r_pb <-  d / sqrt(d^2 + h)
    r_b <- r_pb*(sqrt(p*q)/dnorm(qnorm(p)))
    r <- r_b #r_b = r