How to initialise Data?

[nrgopalrao]

Hi, I have been using SIPPY for system identification as part of my thesis. I have MISO system with discrete time series data of sampling time of 20 ms. I am trying to using either ARMAX or N4SID for this purpose. In both results, the initialisation seems way off the original data but N4SID was comparatively better. For N4SID, after I run the code, and generate the system responses using the results, the initial values seem completely off from the original data. I read in the user manual that "for N4SID’, ’MOESP’ and ’CVA’ methods, x0 is filled with zeros (no estimate of the initial state is made)" but for "fsetSIM.SS_lsim_process_form(A,B,C,D,u,x0)", we can define our own x0. I have made trial and error with this x0 and have been able to fit the data reasonably well but i would like to request you to kindly guide me in understanding how to scientifically choose this x0 or is there a work around during the system identification process itself. The main reason for my doubt is that the original data does not start from 0 but the sys response always try to start from 0 which i think is because the initial states are assumed 0 during system identification. I am a newbie in this field so please excuse my lack of knowledge. I look forward to your response. Thank you Regards, Gopal.

[CPCLAB-UNIPI]

Dear Gopal,

Thanks for your interest in SIPPY. We apologize for the delay in our reply, but we have some organizational problems due to Covid-19 in Italy. As you have read in the user guide no estimate of the initial state is made for methods N4SID ',' MOESP 'and' CVA, but the Parsim methods, for example, ’PARSIM-K’ do it. You may try it.

Where does your data come from? It may be useful to pre-treat your data, for example by scaling your data considering the deviation from a physical equilibrium. ( y_m = y-\bar{y}, u_m=u-\bar{u}) It could be useful for you to read chapter 14 of System Identification: Theory for the User by Lennart Ljung.

Sippy has as optional parameters two centering options: (centering=’MeanVal’ and centering=’InitVal’) You may try to use centering=’InitVal’, and then you add the initial value of your data to the prediction of the identified model.

Best regards, the Sippy teams.

[nrgopalrao]

Hi Sippy teams,

Thank you very much for your response. The PARSIM-K method has worked well for me. In this regard, I just have two more questions.

1. How would i achieve the initialisation for ARMAX method? I have tried to find a solution to it, but unfortunately i have not been successful so far.
2. I was unable to find PARSIM-K as part of Matlab System Identification toolbox. Could you kindly recommend me how to verify these results in Matlab or any other tool that would do the job.

Once again, thank you for responding in such difficult times. I hope you are all safe and healthy.

Best Regards, Gopal

[CPCLAB-UNIPI]

Dear Gopal,

We are happy that Parsim-k has worked well for you. For the ARMAX method, you can try to pre-treat your data as mentioned in the previous answer.

For the verify of the results, you may compare the simulation plot with other State-Space identified model. Try to compare with the N4SID of Matlab System Identification toolbox.

Best Wishes, the Sippy teams.

[nrgopalrao]

Dear SIPPY teams,

I have tried using the ARMAX method on a spring mass system to understand if i am using it correctly. Unfortunately i am not able to understand where exactly I am going wrong. I am unable to reproduce the output from the model created using the input. Input to the system is the force and the output is the position of the spring. I am attaching the script and the excel file for the script. I look forward to your help.

Regards, Gopal springmassdata_200_05302020.xlsx

MIMO_springmass.txt

[CPCLAB-UNIPI]

Dear Gopal,

We apology for the late reply. We hope you solved the problem.

You have some problems with your code. First of all, you have 20000 points but only the first 4000 are different from 0, just use those. (The "GBN" inputs stimulate the system in a better way. You can generate them with the "fset.GBN_seq" function.)

For second you use as "sampling_time = 0.01" but in the simulation you use for the time "npts = len (dataset.index) end_time = npts-1 Time = np.linspace (0, end_time, npts) " In this way you simulate in the system up to a time of "20000" instead of "200".

You should rescale the time "Time = np.linspace (0, end_time, npts) * 0.01"

The last problem is the following: Eout_id1, Time, Xesiml = cnt.lsim (hid1, U_Val [0], Time)

The input of "H" must be "e" and not "u" "Y = G U + H e"

The prediction of the output can also be rewritten in the following form "y_k = H ^ -1 G u_k + (I-H ^ -1) y_k" but the current SIPPY code does not have an integrated function for the conversion . We will probably include it in a future update.

Best Wishes, the Sippy teams.

[nrgopalrao]

Hi SIPPY teams,

Thank you for the response. I have changed a few things like you mentioned and used the scripts. It worked when i used "IC=AIC" etc. It gave me the proper system response but when i use these information criterian, na and nb come out as 2 and 4. But from mathematical description, we know Spring mass damper is a PT2 system so i was expecting na and nb to be 1 and 2. Could you kindly help me understand why we get different values?.

Thank you, Gopal

[CPCLAB-UNIPI]

Dear Gopal,

The identification algorithm does not know the exact order of the system. An information criterion (BIC or AIC) is minimized and the identified system can be of the order of the real system. You can directly choose the order of the identified system if you know the order of the real system, and you want an identified system of that order.

Best Wishes, the Sippy teams.