CODEWITHSUNDEEP
r
Child79
Child79

Reputation: 55

How to write this function to have as output a vector when the input is a vector?

Hello, I want to have in output a vector when the entry is vector without using a function of apply's family. How should I write my function? Thanks. I used this code where I was forced to use two functions.

f1=function(l){
  y= B1 # vector of length N
  li= position # a vector  of length N 
  h=10
  a=(li-l)/h 
  Knorm=dnorm(a)
  b=Knorm*y
  num=sum(b)
  den=sum(Knorm)
  num/den
}
########## Forme vectorielle 
f2 = function(l){
  sapply(l,f1)
} 
L=seq(10000000,11000000,by=1)
f2(L)

If I compute f1(L) I will get one value. That's why I was forced to write a second function to apply my first function to each element of vector L. The purpose is to write it in one function.

Upvotes: 0

Views: 125

Answers (1)

jblood94
jblood94

Reputation: 17001

Use outer and colSums to allow the function to take l as a vector:

f <- function(l){
  y <- B1 # vector of length N
  li <- position # a vector  of length N 
  h <- 10
  a <- outer(li, l, "-")/h
  Knorm <- dnorm(a)
  b <- Knorm*y
  num <- colSums(b)
  den <- colSums(Knorm)
  num/den
}

And here is a simpler equivalent function:

f <- function(l){
  Knorm <- dnorm(outer(position, l, "-")/10)
  colSums(Knorm*B1)/colSums(Knorm)
}

Compare to OP's function:

f1=function(l){
  y= B1 # vector of length N
  li= position # a vector  of length N 
  h=10
  a=(li-l)/h 
  Knorm=dnorm(a)
  b=Knorm*y
  num=sum(b)
  den=sum(Knorm)
  num/den
}

position <- 10:1
B1 <- 1:10

sapply(8:12, f1)
#> [1] 5.300480 5.220937 5.141656 5.062713 4.984177
f(8:12)
#> [1] 5.300480 5.220937 5.141656 5.062713 4.984177

UPDATE

Based on the comments, something like this may work best for the large vectors involved:

library(parallel)

f1 <- function(l) {
  dkAll <- abs(outer(position, l, "-"))
  Knorm <- dnorm(outer(position, l, "-")/pmax(dkAll[order(col(dkAll), dkAll)[seq(70, by = length(position), length.out = length(l))]], 1000))
  colSums(Knorm*y)/colSums(Knorm)
}

y <- seq(1, 100, length.out = 23710)
position <- seq(10351673, 12422082, length.out=23710)
l <- seq(11190000, 11460000, by=10)
# ysmoothed <- f1(l) # memory allocation error
cl <- makeCluster(detectCores())
clusterExport(cl, list("y", "position", "l", 'f1'))
system.time(ysmoothed <- parLapply(cl, l, f1))
#>    user  system elapsed 
#>    0.02    0.00   20.13

Created on 2022-02-02 by the reprex package (v2.0.1)

Upvotes: 1

Related Questions