[R] Regression of complex-valued functions
Rolf Turner
r.turner at auckland.ac.nz
Tue Feb 11 09:42:14 CET 2014
I have not the mental energy to go through your somewhat complicated
example, but I suspect that your problem is simply the following: The
function nls() is trying to minimize a sum of squares, and that does not
make sense in the context of complex observations. That is, nls() is
trying to minimize a sum of terms of the form
(f(x_i,theta) - z_i)^2
where f(x_i,theta) and z_i are complex quantities (and theta is the
vector of parameters with respect to which minimization is "trying to be
done". But as I said, this makes no sense. Not that
(f(x_i,theta) - z_i)^2 is a complex quantity and *not* a positive real
quantity.
You need to minimize a sum of terms of the form
|f(x_i,theta) - z_i|^2
where |w|^2 equal to w times the complex conjugate of w. You get this
in R via the Mod() function, i.e. |w|^2 is (Mod(w))^2.
Thus it seems to me that nls() cannot handle what you want handled. You
need to code up the function to be minimized yourself, using Mod() and
then use optim() or nlm() to minimize it that function with respect to
"theta".
HTH
cheers,
Rolf Turner
On 10/02/14 11:45, Andrea Graziani wrote:
> Hi everyone,
>
> I previously posted this question but my message was not well written and did not contain any code so I will try to do a better job this time.
>
> The goal is to perform a non-linear regression on complex-valued data.
> I will first give a short description of the data and then describe the complex-valued non-linear function.
>
> ***The Data
> Obtained from mechanical tests performed at 5 loading frequencies 0.1, 0.25, 1, 0.4 and 12 Hz, and at 5 temeperatures.
> The independent variable used in the regression is:
> f_data <- rep(c(0.1,0.25,1,4,12),5)
>
> The measured values of the response variable are:
> E <- c(289.7411+ 225.0708i , 386.4417+ 303.5021i , 671.5132+ 521.1253i , 1210.8638+ 846.6253i , 1860.9623+1155.4139i ,
> 862.8984+ 636.2637i , 1159.0436+ 814.5609i , 1919.0369+1186.5679i , 3060.7207+1573.6088i , 4318.1781+1868.4761i ,
> 2760.7782+1418.5450i , 3306.3013+1612.2712i , 4746.6958+1923.8662i , 6468.5769+2148.9502i , 8072.2642+2198.5344i ,
> 6757.7680+2061.3110i , 7591.9787+2123.9168i , 9522.9471+2261.8489i , 11255.0952+2166.6411i , 12601.3970+2120.7178i ,
> 11913.6543+2016.0828i , 12906.8294+2030.0610i , 14343.7693+1893.4877i , 15942.7703+1788.0910i , 16943.2261+1665.9847i)
>
> To visualize the data:
> plot(f_data,Re(E),log="xy")
> plot(f_data,Im(E),log="xy")
> plot(E)
>
> ***Non-linear regression function:
> Obtained from an analytical model
>
> E_2S2P1D <- function(f,logaT,Eg,Ee,k,h,delta,logbeta,logtau)
> Ee+(Eg-Ee)*(
> 1+delta*(2i*pi*10^(logtau)*f*10^logaT)^-k +
> (2i*pi*10^(logtau)*f*10^logaT)^-h +
> (2i*pi*10^(logtau)*f*10^logaT*10^(logbeta))^-1
> )^-1
>
> E_2S2P1D is a complex-valued function (note the imaginary unit "i") where:
> "f" is the real-valued independent variable (i.e. the frequency),
> "logaT","Eg","Ee","k","h","delta","logbeta","logtau" are real-valued parameters.
>
> For example:
>
>> E_2S2P1D(1,0,27000,200,0.16,0.47,1.97,2.2,0.02)
> [1] 9544.759+2204.974i
>
>> E_2S2P1D(1,-1,27000,200,0.16,0.47,1.97,2.2,0.02)
> [1] 6283.748+2088.473i
>
> In order to find the parameters of the regression function I use “nls".
>
> Executing
>
> fit <- nls(E ~ E_2S2P1D(f_data,c(rep(loga40,5),rep(loga30,5),rep(loga20,5),rep(0,5),rep(loga0,5)),
> Eg,Ee,k,h,delta,logbeta,logtau),
> start=list(loga40=-3.76,loga30=-2.63,loga20=-1.39,loga0=1.68,
> Eg=27000,Ee=200,k=0.16,h=0.47,delta=1.97,logbeta=2.2,logtau=-0.02) )
>
> results in a lot of warnings (my translation into english):
> "In numericDeriv(form[[3L]], names(ind), env) :
> imaginary parts removed during the conversion”
>
> I’ve the same error trying:
> y <- E_2S2P1D(f_data,c(rep(-3.76,5),rep(-2.63,5),rep(-1.39,5),rep(0,5),rep(1.68,5)),27000,200,0.16,0.47,1.97,2.2,0.02)
> plot(E,y)
>
> If I got it right R is not able to handle regression problems on complex-valued functions, correct?
> If yes, is there some workaround?
> For example using MSExcel I write the Real and Imaginary parts separately, and than I minimize
>
> error = sum [ ( Re(E) - Re(E_E_2S2P1D) ) / Re(E) )^2 ] + sum [ ( Im(E) - Im(E_E_2S2P1D) ) / Im(E) )^2 ]
>
> using the MSExcel solver function.
>
> thank you for your help
> Andrea
>
> PS
> Thanks to David Winsemius for his suggestions
>
> ______________________________________________
> R-help at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>
More information about the R-help
mailing list