CODEWITHSUNDEEP
rdataframefor-loopode
dingsda345
dingsda345

Reputation: 1

For loops and dataframes, solving differential equation in R

I am having trouble understanding for loops in R. I want to kind of manually solve the equation "exp(kx)=k(x-t)+1" for x. I can't just get x on the left side of the equation alone, so I used ode to get numeric values for exp(k*x). I have estimates for the values of k and t.

library(deSolve)
library(graphics)
#estimated parameter value for k
k<-0.1
#function with ODE system
charnov_intercept=function(t,x,parameters){
  A=parameters[1];
  k=parameters[2];
  traveltime=parameters[3];
  odeVec=rep(0,1);
  odeVec[1]=k*(x[1]-t)+1;
  return(list(odeVec))
}
#give numeric values for exp(k*x) 
numsol_expkT<-lsoda(y=c(0.1),times=seq(0,50,by=0.2),func=charnov_intercept,parms=c(30,0.1,-41))
df<-data.frame(numsol_expkT)

Now I would like to manually solve for x. I want to use a for loop applying ln(exp(k*x), these are the values that should be in df$X1 now) and then dividing this by k to obtain my values for x. I would like to put all of these into new columns in my dataframe.

#ln 
numsol_kT<-vector("double",ncol(df))
for(i in df$X1){
  numsol_kT[i]<-log(i) #achtung, log ist hier nicht log, sondern tatsächlich ln
  print(numsol_kT)
}

This doesn't work, however. R gives me 64 instead of 251 values. I'm brand new at R and have never worked with for loops or anything more than plots before, so I figure I'm doing something fundamentally wrong here. The rest of the code goes:

#put vector as new col in df
df$numsol_kT<-numsol_kT[[2]]

#divide by k and save as new column
for(i in numsol_expkT[,2]){
  numsol_kT[i]<-log(numsol_expkT[i,2])
  print(numsol_kT)
}
df$Tpred_t<-numsol_kT[[2]]

I'm still at the very beginning, but if this doesn't work I can't really continue working on it. Thanks a lot in advance!

Upvotes: 0

Views: 173

Answers (1)

Chris
Chris

Reputation: 2306

Imagining I know nothing about the expected values here and picking up from numsol_expkT:

#give numeric values for exp(k*x) 
numsol_expkT<-lsoda(y=c(0.1),times=seq(0,50,by=0.2),func=charnov_intercept,parms=c(30,0.1,-41))
df<-data.frame(numsol_expkT)
head(df)
  time        X1
1  0.0 0.1000000
2  0.2 0.3020202
3  0.4 0.5040812
4  0.6 0.7061839
5  0.8 0.9083289
6  1.0 1.1105173

ncol(df)
[1] 2
> nrow(df)
[1] 251

it is likely we now want nrow not ncol for numsol_kT

numsol_kT<-vector("double",nrow(df))
for(i in 1:length(df$X1)){
  numsol_kT[i]<-log(df$X1[i]) 
  
} 
str(numsol_kT)
 num [1:251] -2.3026 -1.1973 -0.685 -0.3479 -0.0961 ...
df$numsol_kT <- numsol_kT

head(df)
  time        X1   numsol_kT
1  0.0 0.1000000 -2.30258509
2  0.2 0.3020202 -1.19726137
3  0.4 0.5040812 -0.68501788
4  0.6 0.7061839 -0.34787964

we can make a new receiver in df for next loop, made outside the loop as you did, but in this case a new df column

df$Tpred_t <- NA
for (i in 1:length(numsol_expkT[,2])) {
+ df$Tpred_t[i] <- log(numsol_expkT[i,2])
+ }
> head(df)
  time        X1   numsol_kT     Tpred_t
1  0.0 0.1000000 -2.30258509 -2.30258509
2  0.2 0.3020202 -1.19726137 -1.19726137
3  0.4 0.5040812 -0.68501788 -0.68501788
4  0.6 0.7061839 -0.34787964 -0.34787964
5  0.8 0.9083289 -0.09614875 -0.09614875
6  1.0 1.1105173  0.10482592  0.10482592

And as said, this is just about loop mechanics, not resultant values as I don't know what those should be. for is a counter, i is the count, 1:length( both where are we starting to count and how much. You have to initialize 'receivers' outside the loop, as it is it's own environment, and gone when done, and you did initialize your receivers outside the loop (vector).

Upvotes: 0

Related Questions