How to vectorize integrate when the functions take vectors as arguments?

Question

Hello, I would like to calculate A_i for N=40 individuals. K is the density function of the standard normal distribution. For each of the individuals I have l_i and y_i which are vectors of size 20000. to do this, I first wrote the function that is inside the integral. This function can take vectors as input and return vectors.Then I created a function that calculates the integral. Now suppose that the first 10 individuals have l_i=l1, the next 10 have l_i=l2, the next 10 have l_i=l3 and the last 10 have l_i=l4. And we stack in matrices Y1, Y2, Y3 and Y4 their respective y_i. So Y1, Y2, Y3 and Y4 are matrices of 10 columns.

Now I want to calculate A_i for my 40 individuals.
I used the R code below. But the code takes too long to execute because of the loop. I'm looking for a way to rewrite the code without using a loop. Thank you for helping me.

f1<-function(t,l_i,y_i ){
  Knorm <- dnorm(outer(l_i,t, "-") 
  colSums(Knorm *y_i)/colSums(Knorm)
}
f<-function(t,l_i,y_i){t*f1(t,l_i,y_i)}
Ai_fun<- function(l_i,y_i){
  round(integrate(f,lower = 0,upper = 1,l_i=l_i,y_i=y_i,
                  subdivisions = 2000,stop.on.error = FALSE)$value,4)}
##------------------------------------------------------------------------------ 
##  l1,l2,l3,l4 must be vectors of length=20000
## Y1,Y2,Y3,Y4 must be matrix of 20000 rows and 10 columns
##----------------------------------------------------------------------------------
 
A1<-A2<-A3<-A4<-numeric(10) 
 
for (i  in 1:10) {
 A1[i]<-Ai_fun(l_i = l1,y_i=Y1[,i])  
 A2[i]<-Ai_fun(l_i = l2,y_i=Y2[,i])
 A3[i]<-Ai_fun(l_i = l3,y_i=Y3[,i])
 A4[i]<-Ai_fun(l_i = l4,y_i=Y4[,i])
  
}
A<-c(A1,A2,A3,A4)